Paso 7 (3)- Trayectorias de hospitalización y mortalidad con foco en condiciones vinculadas a trastornos de salud mental y consumo de sustancias posterior a un primer ingreso por alguno de estos trastornos, en usuarios/as jóvenes y adultos emergentes de población general y pertenecientes a pueblos originarios, 2018-2021, Chile

Análisis de sensibilidad para resolución mensual, utilizando aquella solución que obtuvo índices de calidad aceptables.

Configurar

# remover objetos y memoria utilizada
rm(list=ls());gc()

          used (Mb) gc trigger (Mb) max used (Mb)
Ncells  598383 32.0    1303655 69.7   686442 36.7
Vcells 1129253  8.7    8388608 64.0  1876026 14.4

#remover imágenes
while(!dev.cur())dev.off()
cat("\014")

if(Sys.info()["sysname"]=="Windows"){
 folder_path <- ifelse(dir.exists("H:/Mi unidad/PERSONAL ANDRES/UCH_salud_publica/asignaturas/un_inv_II/"),
                       "H:/Mi unidad/PERSONAL ANDRES/UCH_salud_publica/asignaturas/un_inv_II/",
                       ifelse(dir.exists("H:/My drive/PERSONAL ANDRES/UCH_salud_publica/asignaturas/un_inv_II/"),
                       "H:/My drive/PERSONAL ANDRES/UCH_salud_publica/asignaturas/un_inv_II/",
                       "C:/Users/CISS Fondecyt/Mi unidad/Alvacast/SISTRAT 2022 (github)/_proposal_grant/2023/"))
} else {folder_path <- ""}
load(paste0(folder_path,"avance250117_2_25.RData"), verbose = TRUE)

Loading objects:
  tab_clus_reg_pam_om6_q
  surv_diff_4c_neg_sil
  seq_plot_pam_om6_q
  ppoo_clus_pre_pam_om4_q
  suicidio_1_trim_tsm_amplio_2nd
  glosa_pueblo_originario_rec2
  first_episode_tsm_1yr_otro_c6
  episodes_tsm_1yr_tsm_c6
  pamRange_month_om
  second_trim_tus_clas_c6
  second_episode_rec
  suicidio_1_trim_tsm_2nd
  avs_por_cluster_quarter
  dist_quarter_lcs
  km_fit_6c
  lcs_dist_month
  diag_pam_om6_q_6624
  second_com_1_trim_clas_c6
  lcs_dist_quarter
  format_table_vec
  first_episode_com_rec_c6
  tab_clus_sexo_pam_om6_q
  second_episode_com_rec
  tab_clus_prev_pam_om6_q
  ppoo_clus_pre_pam_om4_q_clasificaciones2
  plot_data_melt_pam_om4_q
  second_episode_tsm_1_sem_tsm
  chisq_cramerv
  tab_clas_por_run2
  costmatrix_quarter
  second_sem_tsm_clas
  prim_trim_tsm_clas_c6
  first_episode_rec_c6
  tot_cie10_clas_2nd_c6
  .Random.seed
  episodes_1_sem_tsm
  second_episode_tsm_1_sem_tsm_rec
  episodes_tsm_1yr_otr_c6
  second_episode_tsm_1yr_otro_rec_c6
  best_subset_multinom_interactions_parallel
  tab_clus_compl_pam_om6_q
  second_trim_tsm_1yr_other_clas_c6
  smd_bin
  suicidio_1_trim_tus_amplio_2nd
  ppoo_clus_pre_pam_om4_q_clasificaciones_sin_asw_neg
  tabla_conteo
  second_trim_tsm_clas
  first_episode_rec
  first_episode_1_trim_tsm
  episodes_1_sem_tsm_c6
  new_labels_eng
  in_depresivos
  todes_suicidio_amplio_2nd
  glosa_pueblo_originario_rec_num
  tab_clus_reg_pam_om4_q
  neg_sil_width
  suicidio_1_trim_tsm
  rownames_vec
  dist_month_om
  ppoo_clus_pre_pam_om6_q
  second_trim_com_clas
  sil_neg_pam_om_clus4_q
  ppoo_clus_pre_pam_om4_q_clasificaciones_sin_asw_neg2
  surv_diff_6c
  seq_plot_hc_om2_m
  seqcompare_ppoo_quarter_om
  ppoo_clus_pre_pam_om4_q_clasificaciones
  sil_pam_om_clus6_q_nostd
  categories_pam_om6_q
  first_episode_tsm_1yr_tsm_c6
  pamRange_quarter_lcs
  seqcompare_sex_quarter_om
  todes_suicidio_amplio
  prim_tsm_1_sem_clas_c6
  in_lesiones_amplio
  ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens
  cl
  second_trim_tsm_1yr_tsm_clas_c6
  folder_path
  tabla_combinada
  second_episode_tsm_1yr_tsm_c6
  episodes_1_trim_tsm_c6
  suicidio_1_trim_tsm_amplio
  seqcount_t
  pamRange_quarter_lcs2
  second_episode_tsm_1yr_otro_c6
  second_episode_com_c6
  def_enc17_21
  .QuartoInlineRender
  format_cells
  trim_tasa2_pam_om6_q_cens_cnt
  tabla_prop_fila
  first_episode_tsm_1_sem_tsm_rec_c6
  second_episode_tsm_1_sem_tsm_c6
  first_episode_com_c6
  dist_month_lcs
  ppoo_clus_pre_pam_om4_q_rapa
  surv_obj_4c
  States_Wide.seq_month_t_prim_adm_cens
  tot_cie10_clas_c6
  suicidio_1_sem_tsm_2nd
  rownames_vec_c6
  km_data_4c
  trim_tasa2_pam_om4_q_cens_rate
  tab_cl_mortalidad_pam_om4_q
  prim_sem_tsm_clas
  suicidio_1_sem_tsm_amplio_2nd
  ppoo_clus_pre_pam_om4_q_rapanui
  suicidio_1_trim_tus_2nd
  first_episode_com_rec
  first_episode
  dt_ing_calendar_quarter_t_desde_primera_adm_dedup
  cluster_labels6
  episodes_1_trim_tus
  seqcompare_sex_quarter_lcs
  om_dist_quarter_c
  tab_cluster_region_pam_om4_q
  categories_hac_om2_m
  in_otros_psicoticos
  pamRange_month_lcs2
  sil_neg_pam_om_clus6_q
  States_Wide.seq_month_t_prim_adm_RM_cens
  pamRange_quarter_om2
  cod_reg_homo_pam_om6_q
  first_episode_tsm_1_sem_tsm
  tot_cie10_clas_2nd
  States_Wide.seq_quarter_t_prim_adm_cens
  episodes_1_sem_com_c6
  sil_pam_om_clus4_q_nostd
  first_episode_tsm_1yr_otro_rec_c6
  categories_pam_om4_q
  anova_clus_pam_om4_q
  data_pam_om4_q
  second_episode_1_trim_tsm
  avs_por_cluster_month
  seq_plot_pam_om2_m
  prim_com_1_trim_clas_c6
  seq_plot2_pam_om4_q
  pamRange_month_lcs
  labels_pam_om4_q
  surv_diff_6c_neg_sil
  seqtrate_t
  in_bipolares
  suicidio_1_trim_tus_amplio
  second_episode_tus_rec_c6
  tab_clus_macrozona_pam_om4_q
  table_data_pam_om4_q
  suicidio_com_amplio_2nd
  in_esquizofrenia
  States_Wide.seq_quarter_t_prim_adm_RM_cens
  episodes_1_sem_com
  first_episode_tsm_1_sem_tsm_rec
  costmatrix_month
  ppoo_clus_pre_pam_om4_q_act
  sil_hc_om_clus2_m
  regex_otro
  sil_hc_om_clus2_q
  tab_cl_mortalidad_pam_om6_q
  oneway_anova_effect_size
  States_Wide.seq_quarter_t_prim_adm
  diag_pam_om7_q_6268
  first_episode_tus_c6
  data_melt_pam_om4_q
  prim_trim_tus_clas
  ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens
  lcs_dist_quarter_c
  ratio_plot
  tab_cluster_region_pam_om6_q
  data_long_establecimiento_2024_std
  first_episode_1_trim_tsm_rec_c6
  first_episode_tsm_1yr_tsm_rec_c6
  seq_plot_pam_om4_q_eng
  seq_plot_pam_om2_q
  dt_ing_calendar_month_t_desde_primera_adm_dedup
  cod_reg_homo_pam_om4_q
  om_dist_month_c
  todes_suicidio
  episodes_1_trim_tsm
  surv_diff_4c
  second_episode_tsm_1_sem_tsm_rec_c6
  anova_clus_pam_om6_q
  States_Wide.seq_month_t_prim_adm
  suicidio_1_sem_tsm_amplio
  second_episode_com_rec_c6
  km_fit_4c
  second_episode_tsm_1yr_tsm_rec_c6
  df_filled2
  om_dist_month
  trim_tasa_pam_om6_q_cens_rate
  seq_plot2_pam_om6_q
  lcs_dist_month_c
  agregar_columnas_icd
  om_dist_quarter
  sil_pam_om_clus4_q
  first_episode_1_trim_tsm_c6
  seqcompare_ppoo_quarter_lcs
  new_labels
  glosa_pueblo_originario_rec
  sil_pam_om_clus2_q
  tab_clus_macrozona_pam_om6_q
  dist_quarter_om
  seqcompare_ppoo_quarter_om_ppoo_rec
  suicidio_1_sem_tsm
  sens_min_anios_clus_pam_om6_q
  new_labels2
  seq_plot_pam_om3_q
  seq_plot_hc_om3_m
  prim_trim_tsm_1yr_tsm_clas_c6
  prim_trim_com_clas
  episodes_1_trim_tus_c6
  max_clusters
  new_labels_eng2
  second_episode_tus
  best_subset_multinom
  second_episode_1_trim_tsm_c6
  suicidio_com
  pamRange_month_om2
  first_episode_tsm_1_sem_tsm_c6
  sil_pam_om_clus3_q
  best_subset_multinom_interactions
  tab_clas_por_run
  diag_pam_om4_6522
  as.data.frame.TableOne
  second_episode_tus_c6
  first_episode_tus
  in_substance
  resultados_list
  episodes_1_overall
  sil_pam_om_clus2_m
  sil_hc_om_clus3_m
  seq_plot_pam_om4_q
  in_personalidad
  ppoo_clus_pre_pam_om6_q_rapanui
  prim_trim_tsm_1yr_other_clas_c6
  prim_trim_tsm_clas
  prim_trim_tus_clas_c6
  first_episode_com
  in_lesiones
  .main
  second_episode_1_trim_tsm_rec
  ppoo_inclusivo
  multinom_pivot_wider
  todes_suicidio_2nd
  num_cores
  second_tsm_1_sem_clas_c6
  sil_pam_om_clus6_q
  cluster_labels
  pamRange_quarter_om
  second_trim_tus_clas
  second_trim_tsm_clas_c6
  seq_mean_t
  tab_clus_sexo_pam_om4_q
  filt_data_nueva
  in_organic
  tab_clus_prev_pam_om4_q
  suicidio_com_2nd
  diag_pam_om7_q_6522
  fisher_posthoc_all_cols
  sens_min_anios_clus_pam_om4_q
  tabla_prop_col
  trim_tasa2_pam_om6_q_cens_rate
  second_episode_tus_rec
  surv_obj_6c
  second_episode
  trim_tasa_pam_om4_q_cens_cnt
  second_episode_rec_c6
  in_ansiedad
  first_episode_tus_rec_c6
  km_data_6c
  suicidio_com_amplio
  suicidio_1_trim_tus
  second_episode_com
  trim_tasa2_pam_om4_q_cens_cnt
  first_episode_1_trim_tsm_rec
  pattern
  tab_clus_compl_pam_om4_q
  mort_rate_tus_resto
  trim_tasa_pam_om6_q_cens_cnt
  trim_tasa_pam_om4_q_cens_rate
  tot_cie10_clas
  first_episode_tus_rec
  second_episode_1_trim_tsm_rec_c6

Paquetes estadísticos

#elegir repositorio
if(Sys.info()["sysname"]=="Windows"){
  options(repos = c(CRAN = "https://cran.dcc.uchile.cl/"))
}
options(install.packages.check.source = "yes") # Chequea la fuente de los paquetes

#borrar caché
#system("fc-cache -f -v")

#elegir repositorio
if(Sys.info()["sysname"]=="Windows"){
  options(repos = c(CRAN = "https://cran.dcc.uchile.cl/"))
}
options(install.packages.check.source = "yes") # Chequea la fuente de los paquetes

#borrar caché
#system("fc-cache -f -v")

if(!require(pacman)){install.packages("pacman");require(pacman)}

pacman::p_unlock(lib.loc = .libPaths()) #para no tener problemas reinstalando paquetes

if(Sys.info()["sysname"]=="Windows"){
if (getRversion() != "4.4.1") { stop("Requiere versión de R 4.4.1. Actual: ", getRversion()) }
}

cat("quarto version: "); system("quarto --version") 

quarto version: 

[1] 0

if(!require(job)){install.packages("job");require(job)}
if(!require(kableExtra)){install.packages("kableExtra");require(kableExtra)}
if(!require(tidyverse)){install.packages("tidyverse");require(tidyverse)}
if(!require(cluster)){install.packages("cluster"); require(cluster)}
if(!require(WeightedCluster)){install.packages("WeightedCluster"); require(WeightedCluster)}
if(!require(devtools)){install.packages("devtools"); require(devtools)}
if(!require(TraMineR)){install.packages("TraMineR"); require(TraMineR)}
if(!require(TraMineRextras)){install.packages("TraMineRextras"); require(TraMineRextras)}
if(!require(NbClust)){install.packages("NbClust"); require(NbClust)}
if(!require(haven)){install.packages("haven"); require(haven)}
if(!require(ggseqplot)){install.packages("ggseqplot"); require(ggseqplot)}
if(!require(grid)){install.packages("grid"); require(grid)}
if(!require(gridExtra)){install.packages("gridExtra"); require(gridExtra)}
if(!require(Tmisc)){install.packages("Tmisc"); require(Tmisc)}
if(!require(factoextra)){install.packages("factoextra"); require(factoextra)}
if(!require(stargazer)){install.packages("stargazer"); require(stargazer)}
if(!require(gtsummary)){install.packages("gtsummary"); require(gtsummary)}
if(!require(lmtest)){install.packages("lmtest"); require(lmtest)}
if(!require(emmeans)){install.packages("emmeans"); require(emmeans)}
if(!require(fpp2)){install.packages("fpp2"); require(fpp2)}
if(!require(purrr)){install.packages("purrr"); require(purrr)}
if(!require(forecast)){install.packages("forecast"); require(forecast)}
if(!require(magrittr)){install.packages("magrittr"); require(magrittr)}
if(!require(foreach)){install.packages("foreach"); require(foreach)}
if(!require(doParallel)){install.packages("doParallel"); require(doParallel)}
if(!require(progressr)){install.packages("progressr"); require(progressr)}
if(!require(chisq.posthoc.test)){devtools::install_github("ebbertd/chisq.posthoc.test")}
if(!require(rstatix)){install.packages("rstatix"); require(rstatix)}
if(!require(rio)){install.packages("rio"); require(rio)}
if(!require(cowplot)){install.packages("cowplot"); require(cowplot)}
if(!require(DiagrammeR)){install.packages("DiagrammeR"); require(DiagrammeR)}
if(!require(DiagrammeRsvg)){install.packages("DiagrammeRsvg"); require(DiagrammeRsvg)}
if(!require(rsvg)){install.packages("rsvg"); require(rsvg)}
if(!require(survminer)){install.packages("survminer"); require(survminer)}
if(!require(epitools)){install.packages("epitools"); require(epitools)}
if(!require(biostat3)){install.packages("biostat3"); require(biostat3)}

Resultados

2. Mensual

2.1. Sensibilidad= PAM (OM), sol 2 cluster- diagnósticos

invisible("Me da buena: 0,61 en promedio. Se mantiene. El problema está con 5710 es negativo")
sil_pam_om_clus2_m_nostd<-
  silhouette(as.integer(pamRange_month_om$clustering$cluster2), as.dist(dist_month_om))

# Crear etiquetas personalizadas
cluster_labels2_m <- paste0("Cluster ", 
          seq_along(attr(summary(sil_pam_om_clus2_m_nostd)$clus.avg.widths, "dimnames")[[1]]), 
          ":\nAWS ", sprintf("%1.2f",summary(sil_pam_om_clus2_m_nostd)$clus.avg.widths))

# Graficar con etiquetas personalizadas
fviz_silhouette(
  sil_pam_om_clus2_m_nostd, 
  lab.clusters = cluster_labels2_m, # Etiquetas personalizadas para los clústeres
  print.summary=F) +
  scale_fill_grey(start = 0.2, end = 0.8, labels = cluster_labels2_m) +  # Escala de grises
  scale_color_grey(start = 0.2, end = 0.8, labels = cluster_labels2_m)+   # Escala de grises para los bordea
  ggtitle(NULL)+
  labs(y="Ancho medio de la silueta", x="Conglomerados")# Elimina el título

ggsave("_figs/sil_plot_pam_om2_m_25.png", width = 8, height = 5, dpi = 500)

# Crear la tabla de frecuencias proporcionales redondeada
tabla_proporciones2 <- round(prop.table(table(pamRange_month_om$clustering$cluster2, 
                                             pamRange_quarter_om$clustering$cluster4), 2), 2)

# Convertir la tabla a un formato limpio con kable
knitr::kable(tabla_proporciones2, 
             caption = "Proporciones de Clusters, Solución de 2 vs. 4 conglomerados", 
             col.names = c("5939, Un semestre TSM(1)", "5989, Comorbilidad un trimestre(2)", 
                           "6025, Un trimestre, TUS(3)", "6035, Un trimestre, TSM(4)"), 
             align = "c")

Proporciones de Clusters, Solución de 2 vs. 4 conglomerados

	5939, Un semestre TSM(1)	5989, Comorbilidad un trimestre(2)	6025, Un trimestre, TUS(3)	6035, Un trimestre, TSM(4)
6612	0.01	0.07	1	0
6623	0.99	0.93	0	1

tabla_conteo2 <- table(pamRange_quarter_om$clustering$cluster2, 
                      pamRange_quarter_om$clustering$cluster4)

# Proporción por fila
tabla_prop_fila2 <- prop.table(tabla_conteo2, 1)

# Proporción por columna
tabla_prop_col2 <- prop.table(tabla_conteo2, 2)

# Redondear ambas

tabla_prop_fila2 <- sprintf("%1.1f", tabla_prop_fila2*100)
tabla_prop_col2 <- sprintf("%1.1f", tabla_prop_col2*100)


# Combinar ambas proporciones en el formato "fila / col"
tabla_combinada2 <- matrix(
  paste0(tabla_prop_fila2, " / ", tabla_prop_col2),
  nrow = nrow(tabla_conteo2),
  ncol = ncol(tabla_conteo2),
  dimnames = dimnames(tabla_conteo2)
)


# Convertir la tabla a un formato limpio con kable
tabla_combinada2 |> 
  kable(
    caption = "Proporciones de Clusters (por fila / por columna), Solución de 6 vs. 4 conglomerados", 
    col.names = c("6522, Un semestre TSM(2)", 
                  "6574, Comorbilidad un trimestre(3)", 
                  "6612, Un trimestre, TUS(4)", 
                  "6623, Un trimestre, TSM(4)"), 
    align = "c"
  ) |> 
  kable_classic() |> 
  footnote(
    general = "6612= Primer mes, TUS; 6623= Primer mes, TSM",
    general_title = "Nota:",
    footnote_as_chunk = TRUE,
    escape = FALSE
  )

Proporciones de Clusters (por fila / por columna), Solución de 6 vs. 4 conglomerados

	6522, Un semestre TSM(2)	6574, Comorbilidad un trimestre(3)	6612, Un trimestre, TUS(4)	6623, Un trimestre, TSM(4)
6612	0.3 / 0.5	0.8 / 2.6	98.9 / 100.0	0.0 / 0.0
6623	6.8 / 99.5	3.8 / 97.4	0.0 / 0.0	89.4 / 100.0
Nota: 6612= Primer mes, TUS; 6623= Primer mes, TSM

cat("2025-03-01: general")

2025-03-01: general

first_episode<-
df_filled2 |> 
    dplyr::filter(run %in% ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$run)|> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 1)|>
    dplyr::ungroup()

second_episode<-
df_filled2 |> 
    dplyr::filter(run %in% ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$run)|> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 2)|>
    dplyr::ungroup()


first_episode_rec_c2<-
agregar_columnas_icd(first_episode)
  
tot_cie10_clas_c2 <-
rbind.data.frame(
  format_table_vec(table(factor(first_episode_rec_c2$organic, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$substance, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$esquizofrenia, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$otros_psicoticos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$bipolares, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$depresivos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$ansiedad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$personalidad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$lesiones, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_rec_c2$lesiones_amplio, levels = c(FALSE, TRUE))))
  )|>
    (\(df) {
  colnames(df) <- c("FALSE", "TRUE"); rownames(df)<-NULL; df
      })()
#         FALSE        TRUE
# 1 4858 (80.5) 1180 (19.5)
# 2 5612 (92.9)   426 (7.1)
# 3 5645 (93.5)   393 (6.5)
# 4 5652 (93.6)   386 (6.4)
# 5 4215 (69.8) 1823 (30.2)
# 6 5072 (84.0)  966 (16.0)
# 7 4903 (81.2) 1135 (18.8)
# 

second_episode_rec_c2<-
agregar_columnas_icd(second_episode)
  
tot_cie10_clas_2nd_c2<-
rbind.data.frame(
  format_table_vec(table(factor(second_episode_rec_c2$organic, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$substance, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$esquizofrenia, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$otros_psicoticos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$bipolares, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$depresivos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$ansiedad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$personalidad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$lesiones, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_rec_c2$lesiones_amplio, levels = c(FALSE, TRUE))))
  )|>
    (\(df) {
  colnames(df) <- c("FALSE", "TRUE"); rownames(df)<-NULL; df
      })()

cat(paste0("Second episodes: ", 
df_filled2 |> 
    dplyr::filter(run %in% ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 2)|>
    dplyr::ungroup() |> nrow()
),"\n")

Second episodes: 3117 

#:#:#:#:#:#:#:#:#:#:
cat("2025-04-26: TSM 1 semestre")

2025-04-26: TSM 1 semestre

first_episode_1_month_tsm_c2<-
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6623", clus_pam_om2))$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 1)|>
    dplyr::ungroup()
second_episode_1_month_tsm_c2<-
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6623", clus_pam_om2))$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 2)|>
    dplyr::ungroup()

cat(paste0("Second episodes, First quarter, MH: ", nrow(second_episode_1_month_tsm_c2),"\n"))

Second episodes, First quarter, MH: 2805

#:#:#:#:#:#:#:#:#:#:
cat("2025-04-28: TUS")

2025-04-28: TUS

first_episode_tus_c2<-
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6612", clus_pam_om2))$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 1)|>
    dplyr::ungroup()
second_episode_tus_c2<-
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6612", clus_pam_om2))$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::filter(row_number() == 2)|>
    dplyr::ungroup()

cat(paste0("Second episodes, First quarter, SUD: ", nrow(second_episode_tus_c2),"\n"))

Second episodes, First quarter, SUD: 312

#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:

first_episode_1_mnth_tsm_rec_c2<-
agregar_columnas_icd(first_episode_1_month_tsm_c2)

prim_mnth_tsm_clas_c2<-
rbind.data.frame(
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$organic, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$substance, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$esquizofrenia, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$bipolares, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$depresivos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$ansiedad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$otros_psicoticos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$personalidad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$lesiones, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_1_mnth_tsm_rec_c2$lesiones_amplio, levels = c(FALSE, TRUE))))
)|>
    (\(df) {
  colnames(df) <- c("FALSE", "TRUE"); rownames(df)<-NULL; df
      })()

second_episode_1_mnth_tsm_rec_c2<-
agregar_columnas_icd(second_episode_1_month_tsm_c2)

second_mnth_tsm_clas_c2 <-
rbind.data.frame(
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$organic, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$substance, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$esquizofrenia, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$otros_psicoticos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$bipolares, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$depresivos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$ansiedad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$personalidad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$lesiones, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_1_mnth_tsm_rec_c2$lesiones_amplio, levels = c(FALSE, TRUE))))
)|>
    (\(df) {
  colnames(df) <- c("FALSE", "TRUE"); rownames(df)<-NULL; df
      })()

#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:

first_episode_tus_rec_c2<-
agregar_columnas_icd(first_episode_tus_c2)#first_episode_com_c6

prim_trim_tus_clas_c2<-
rbind.data.frame(
  format_table_vec(table(factor(first_episode_tus_rec_c2$organic, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$substance, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$esquizofrenia, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$otros_psicoticos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$bipolares, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$depresivos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$ansiedad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$personalidad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$lesiones, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(first_episode_tus_rec_c2$lesiones_amplio, levels = c(FALSE, TRUE))))
)|>
    (\(df) {
  colnames(df) <- c("FALSE", "TRUE"); rownames(df)<-NULL; df
      })()

second_episode_tus_rec_c2<-
agregar_columnas_icd(second_episode_tus_c2)#

second_trim_tus_clas_c2<-
rbind.data.frame(
  format_table_vec(table(factor(second_episode_tus_rec_c2$organic, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$substance, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$esquizofrenia, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$otros_psicoticos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$bipolares, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$depresivos, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$ansiedad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$personalidad, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$lesiones, levels = c(FALSE, TRUE)))),
  format_table_vec(table(factor(second_episode_tus_rec_c2$lesiones_amplio, levels = c(FALSE, TRUE))))
)|>
    (\(df) {
  colnames(df) <- c("FALSE", "TRUE"); rownames(df)<-NULL; df
      })()

#::::::::##::::::::##::::::::##::::::::##::::::::##::::::::##::::::::#
#::::::::##::::::::##::::::::##::::::::##::::::::##::::::::##::::::::##::::::::#
#::::::::##::::::::##::::::::##::::::::##::::::::##::::::::##::::::::#
#
#

cbind.data.frame(Total= tot_cie10_clas_c2[,2],
                 Total_2nd= tot_cie10_clas_2nd_c2[,2],
                 Prim_trim_tsm= prim_mnth_tsm_clas_c2[,2],
                 Prim_trim_tsm_2nd= second_mnth_tsm_clas_c2[,2],
                 Prim_trim_tus= prim_trim_tus_clas_c2[,2],
                 Prim_trim_tus_2nd= second_trim_tus_clas_c2[,2]
           )|>
    (\(df) {
  rownames(df) <- rownames_vec_c6
  rio::export(df,"_output/diagnostic_classification_cie10_c2.xlsx")
  df
      })() |> 
  knitr::kable("markdown", caption= "Clasificación diagnóstica CIE-10")

Clasificación diagnóstica CIE-10

	Total	Total_2nd	Prim_trim_tsm	Prim_trim_tsm_2nd	Prim_trim_tus	Prim_trim_tus_2nd
Trastornos mentales orgánicos: F00.0-F09.9	100 (1.5)	19 (0.6)	100 (1.7)	17 (0.6)	0 (0)	2 (0.6)
Trastornos por uso de sustancias: F10.0-F19.9	1293 (19.5)	368 (11.8)	549 (9.4)	258 (9.2)	744 (96.7)	110 (35.3)
Esquizofrenia: F20.0-F20.9	441 (6.7)	216 (6.9)	441 (7.5)	211 (7.5)	0 (0)	5 (1.6)
Otros trastornos psicóticos no afectivos: F21.0-F29.9	423 (6.4)	115 (3.7)	411 (7.0)	106 (3.8)	1 (0.1)	9 (2.9)
Trastornos bipolares: F30.0-F31.9	411 (6.2)	199 (6.4)	2007 (34.3)	195 (7.0)	0 (0)	4 (1.3)
Trastornos depresivos y otros del estado de ánimo: F32.0-F39.9	2009 (30.3)	470 (15.1)	1082 (18.5)	463 (16.5)	2 (0.3)	7 (2.2)
Trastornos de ansiedad: F40.0-F49.9	1082 (16.3)	182 (5.8)	422 (7.2)	173 (6.2)	0 (0)	9 (2.9)
Trastornos de la personalidad: F60.0-F69.9	1234 (18.6)	334 (10.7)	1231 (21.0)	320 (11.4)	3 (0.4)	14 (4.5)
Heridas autoinflingidas (info. intrahosp): X60-X84	237 (3.6)	152 (4.9)	220 (3.8)	148 (5.3)	17 (2.2)	4 (1.3)
Intoxicación/efectos tóxicos de químicos y sustancias: T36-T65	288 (4.3)	170 (5.5)	253 (4.3)	164 (5.8)	35 (4.6)	6 (1.9)

episodes_1_mnth_overall<-
df_filled2 |> 
  dplyr::filter(run %in% ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::summarise(n= n())|>
    dplyr::ungroup()


episodes_1_mnth_tsm_c2<-
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6623", clus_pam_om2))$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::summarise(n= n())|>
    dplyr::ungroup()

episodes_1_mnth_tus_c2<-
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6612", clus_pam_om2))$run) |> 
    dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) |> 
    dplyr::group_by(run) |>
    dplyr::summarise(n= n())|>
    dplyr::ungroup()

#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
rbind.data.frame(
df_filled2 |> 
  dplyr::filter(run %in% ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$run) |> 
  pull(days_elapsed) |> psych::describe(quant=c(.25,.75)),
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6623", clus_pam_om2))$run) |> 
  pull(days_elapsed) |> psych::describe(quant=c(.25,.75)),
df_filled2 |> 
  dplyr::filter(run %in% subset(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, grepl("6612", clus_pam_om2))$run) |> 
  pull(days_elapsed) |> psych::describe(quant=c(.25,.75))
)|>
  (\(df) {
 rownames(df) <- c("Total", "M1, TSM", "M1, TUS"); df
      })() |> 
  dplyr::select(-any_of(c("trimmed", "mad", "range", "skew", "kurtosis","se"))) |> 
      (\(df) {
  if (interactive()) {df|> rio::export("clipboard")}
  knitr::kable(df, "markdown", caption="Días de hospitalización por cluster", digits=2)
  })()

Días de hospitalización por cluster

	vars	n	mean	sd	median	min	max	Q0.25	Q0.75
Total	1	13230	13.41	30.29	6	0	1425	2	16
M1, TSM	1	11747	13.78	31.08	6	0	1425	2	17
M1, TUS	1	1483	10.47	22.95	3	0	365	1	11

rbind.data.frame(
psych::describe(episodes_1_mnth_overall$n, quant=c(.25,.75)),
psych::describe(episodes_1_mnth_tsm_c2$n, quant=c(.25,.75)),
psych::describe(episodes_1_mnth_tus_c2$n, quant=c(.25,.75))
) |> 
  (\(df) {
 rownames(df) <- c("Total", "M1, TSM", "M1, TUS")
   dplyr::select(df, -any_of(c("trimmed", "mad", "range", "skew", "kurtosis","se")))
  })() |> 
  (\(df) {
  if (interactive()) {df|> rio::export("clipboard")}
  knitr::kable(df, "markdown", caption="Número de episodios por cluster", digits=2)
  })()

Número de episodios por cluster

	vars	n	mean	sd	median	min	max	Q0.25	Q0.75
Total	1	6626	2.00	1.75	1	1	36	1	2
M1, TSM	1	5857	2.01	1.73	1	1	36	1	2
M1, TUS	1	769	1.93	1.94	1	1	21	1	2

janitor::chisq.test(df_filled2[,c("run","glosa_pueblo_originario")]|> 
        dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2")], by="run", multiple="first")|> 
        janitor::tabyl(glosa_pueblo_originario, clus_pam_om2))

    Pearson's Chi-squared test

data:  janitor::tabyl(dplyr::left_join(df_filled2[, c("run", "glosa_pueblo_originario")],     ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,         c("run", "clus_pam_om2")], by = "run", multiple = "first"),     glosa_pueblo_originario, clus_pam_om2)
X-squared = 19.966, df = 8, p-value = 0.01046

# X-squared = 19.966, df = 8, p-value = 0.01046

Generamos un gráfico de PPOO por cada conglomerado.

# 2025-05-03: cambiar etiquetas de los dos grupos
levels(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2) <- sub("Un trimestre", "Primer mes", levels(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2))
levels(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2) <- sub("TSM\\(1\\)", "TUS(1)", levels(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2))

#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_


ppoo_clus_pre_pam_om2_m<-
  df_filled2[,c("run","glosa_pueblo_originario")]|> 
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2")], by="run", multiple="first")|> 
  dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "ci_cya_conadi", "ppoo_minsal_y_rsh_2010")], by="run", multiple="first")|> 
  dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario == "NINGUNO" & (ci_cya_conadi ==1 | ppoo_minsal_y_rsh_2010==1)~ "DESCONOCIDO", T~glosa_pueblo_originario)) |> 
  janitor::tabyl(glosa_pueblo_originario_rec, clus_pam_om2) |> 
  janitor::adorn_percentages("row")

#0 pacientes en cluster TUS YAGÁN (YÁMANA), AYMARA, COLLA, DIAGUITA, KAWÉSQAR; 
#Pocos pacientes en cluster TUS: OTRO=2, RAPA NUI= 5
#Pocos pacientes en TSM: YAGAN=2, DIAGUITA=3, KAWESQAR=4


invisible("2025-05-03")
cat("Sólo RAPA NUI y MAPUCHE (ver por RUN, al menos una observación que se identifique RAPA NUI o MAPUCHE)\n")
ppoo_clus_pre_pam_om2_m_act<-
  df_filled2[,c("run","glosa_pueblo_originario")]|> 
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2")], by="run", multiple="first")|> 
  dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "ci_cya_conadi", "ppoo_minsal_y_rsh_2010")], by="run", multiple="first")|> 
  dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario == "NINGUNO" & (ci_cya_conadi ==1 | ppoo_minsal_y_rsh_2010==1)~ "DESCONOCIDO", T~glosa_pueblo_originario)) |> 
  group_by(run, clus_pam_om2)|>
  summarise(sum_ppoo= sum(grepl("RAPA|MAPU", glosa_pueblo_originario_rec)), .groups="drop")|>
  ungroup()|> 
  mutate(sum_bin_ppo= ifelse(sum_ppoo>0, 1, 0))|>
  janitor::tabyl(sum_bin_ppo, clus_pam_om2)
janitor::chisq.test((ppoo_clus_pre_pam_om2_m_act))
#X-squared = 4.2837, df = 1, p-value = 0.03848 # ES SIGNIFICATIVO!, x2 OR de pertenecer a uno de estos pueblos en TUS; 1.5 veces la prevalencia

cat("Sólo RAPA NUI (ver por RUN, al menos una observación que se identifique RAPA NUI)\n") #NO SIGNIFICATIVO
ppoo_clus_pre_pam_om2_m_rapa<-
  df_filled2[,c("run","glosa_pueblo_originario")]|> 
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2")], by="run", multiple="first")|> 
  dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "ci_cya_conadi", "ppoo_minsal_y_rsh_2010")], by="run", multiple="first")|> 
  dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario == "NINGUNO" & (ci_cya_conadi ==1 | ppoo_minsal_y_rsh_2010==1)~ "DESCONOCIDO", T~glosa_pueblo_originario))|> 
  group_by(run, clus_pam_om2)|>
  summarise(sum_ppoo= sum(grepl("RAPA", glosa_pueblo_originario_rec)), .groups="drop")|>
  ungroup()|> 
  mutate(sum_bin_ppo= ifelse(sum_ppoo>0, 1, 0))|>
  janitor::tabyl(sum_bin_ppo, clus_pam_om2)

janitor::fisher.test((ppoo_clus_pre_pam_om2_m_rapa))
#p-value = 0.3708
#OR:   1.694165 
cat("son 11 casos!\n")


cat("Clasificando por Mapuche, otro no-mapuche, ninguno puro y desconocido puro")
ppoo_clus_pre_pam_om2_m_clasificaciones<-
df_filled2[, c("run", "glosa_pueblo_originario")]|>
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2")], by="run", multiple="first")|> 
  dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "ci_cya_conadi", "ppoo_minsal_y_rsh_2010")], by="run", multiple="first")|>
  mutate(glosa_pueblo_originario_rec = if_else(
    glosa_pueblo_originario == "NINGUNO" & (ci_cya_conadi ==1 | ppoo_minsal_y_rsh_2010==1), "DESCONOCIDO",
    glosa_pueblo_originario
  ))|>
  group_by(run, clus_pam_om2)|>
  summarise(
    tiene_otro_po = any(grepl("AYMARA|COLLA|DIAGUITA|KAWÉSQAR|RAPA NUI|YAGÁN|OTRO", glosa_pueblo_originario_rec)),
    tiene_mapuche = any(grepl("MAPUCHE", glosa_pueblo_originario_rec)),
    tiene_ninguno = any(glosa_pueblo_originario_rec == "NINGUNO"),
    tiene_desconocido = any(glosa_pueblo_originario_rec == "DESCONOCIDO"),
    .groups = "drop"
  ) |>
  mutate(clas_por_run = case_when(
    tiene_mapuche ~ "Mapuche",
    !tiene_mapuche & tiene_otro_po ~ "Otros reportados (no mapuche)",
    tiene_ninguno & !tiene_otro_po & !tiene_mapuche ~ "Ninguno puro",
    tiene_desconocido & !tiene_ninguno & !tiene_otro_po & !tiene_mapuche ~ "Desconocido puro", T~NA_character_
  )) |>
  janitor::tabyl(clas_por_run, clus_pam_om2, show_na = T)

janitor::adorn_percentages(ppoo_clus_pre_pam_om4_q_clasificaciones, "col")|> 
    mutate_if(is.numeric, ~round(.*100, 1)) |> mutate(clas_por_run= factor(clas_por_run, levels= c("Ninguno puro", "Desconocido puro", "Otros reportados (no mapuche)", "Mapuche")))|> 
  arrange(clas_por_run)|>
  (\(df) {
  if (interactive()) {df|> rio::export("clipboard")}
  knitr::kable(df, caption="Porcentajes por columna, conglomerado vs. PPOO")
  })()

janitor::fisher.test(ppoo_clus_pre_pam_om2_m_clasificaciones, simulate.p.value=T, B=1e5)
#p-value = 0.001


chisq.posthoc.test(ppoo_clus_pre_pam_om2_m_clasificaciones[,2:3],simulate.p.value=T, B=1e5)|> 
  (\(df) { dplyr::mutate_at(df, 3:ncol(df), ~ round(as.numeric(gsub("\\*", "", .)), 3))
  })()|>
  group_by(Dimension)|>
  (\(df) { summarise(df, across(2:(ncol(df)-1), ~ paste0(first(sprintf("%1.2f", .)), " (p=", last(sprintf("%1.3f", .)), ")")))
  })()|>
  (\(df) { dplyr::mutate_at(df, 2:length(names(df)), ~ gsub("p\\=0.000)", "p<0.001)", .))
  })()|>
  mutate(Dimension= ppoo_clus_pre_pam_om4_q_clasificaciones[,1]) |> 
  mutate(Dimension= factor(Dimension, levels= c("Ninguno puro", "Desconocido puro", "Otros reportados (no mapuche)", "Mapuche"))) |> 
  arrange(Dimension) |>
  mutate(Dimension= as.character(Dimension)) |> 
  knitr::kable("markdown", caption = "Comparación post-hoc, conglomerado-PPOO autorreporte")

cat("Tampoco resulta significativo en https://adaptivedesignstrial.shinyapps.io/posthoc/\n")


invisible("Otros reportados (no mapuche) pertenecen menos Primer mes TUS(1); The Bonferroni corrected significant level for Exact approach is: 0.05/(R*C)=
0.00625; pero no habia nada destacado en verde o amarillo")

invisible("pvalue based on Adjusted residual y pvalue based on exact approach, son significativos al 0.01 par Otros reportados y baja pertenencia  TUS; pero no hay yellow or green are significant based on the exact approach under the Holm-Bonferroni correction.")

invisible("0.04 tiene los pvalues basados en residuos ajustados y pvalues basados en la aproximación exacta , en el que Mapuches tienen más proporción TUS")
invisible("0.05 tiene los pvalues basados en residuos ajustados y pvalues basados en la aproximación exacta , en el que Desconocido puro tienen más proporción TUS")

invisible("Los desconocido puro son quienes tienen i_cya_conadi ==1 | ppoo_minsal_y_rsh_2010==1, es decir, históricamente son reconocidos o se han autoidentificado, pero no en MINSAL (que es el que registra la mierdita)")

#####::::::::::::::::::::::::::::::::::::::::::::::::::::::
#####::::::::::::::::::::::::::::::::::::::::::::::::::::::
#####::::::::::::::::::::::::::::::::::::::::::::::::::::::
cat("Ahora descartando ASWs negativos")

run_aws_neg_m_c2<- 
c("01a17af6c7b8aaba09226eff980bda2c6257c7a830e6be8e853a67743070c779", 
"05f17cc6fc0a8e7e45687470ce81ac223f92c6c4320e2cc293d3fa853d1d4478", 
"06129e059d0ef701be2ae0ff8b4c1f604fe2286efc1b05db190b7844a3cb4879", 
"0c1ed1bb2ff17afe1bd9ccc41f3a3028750a6549a056ab657a8be5d231ca1452", 
"0c8f79a7c405104ce14e8d7b31dcc6fb70248fd6ae7de55a642d30eb12d86721", 
"10f8b961d7cbff7cec8e4396f7a3936e54251b78c82c6ce15c84398ded89d845", 
"143139a49434bcf76fdc0cb3d94decc1070e4d229351c0a49e19bb4dc5ae4d6e", 
"156fb1388cc49601df8d578e6d2c04ef52a9908cb5dc5c1392d4a726223c1fa3", 
"15e3dd4ec5b6fe40be3c8ca7b1bfc648c502d6a4488bd65091e32c751c747872", 
"15ebcf907feeff230ced69bc57cff0b53c02af27aeb552fedf80e40e7d3b8639", 
"16d0fda9b29b8afa26be40a49accfd7ecc0eb8f4e41d7b0ec4e5cb791ca23622", 
"19e0e606397d62e5b9d4d71797877302588c7e17c46b1ed4b9712f11910b9677", 
"1c45c5a100a9bf4bf31d48c1794295388f766a1d6916c69d0e8bf748843d37ad", 
"1caba8fbe7076a293a2ae24160ddd627b6a89aea6271cc13f1b9c4f9f48b1b54", 
"1dd360fb4101d4f4b78162344d80b257c745d2eb81d57f96526a29948b31864c", 
"1f786ae1dcea0aebc6339e03176f644a248660bae766fa651f4d4f6bf0eb3c0e", 
"2059f8a5b4ff2fca074ec1e0b549dcccde4231e77a51943b4bac770ea1c70439", 
"21d49f5e2ead538e7802f5854106dc7da7f200f538d13100bad8e9be9da99317", 
"22d8751072f17d01fed159cf7b26d0128fcae823732b484a3409593e508a0321", 
"22d8a99fb68720cceb1f5c342017e3298f37c2e3d120b0f7a63d67216b0d501c", 
"22f44ce402659a7c087bf0b04b38480e8ec81f727905cae84db4897ffa5979c8", 
"230c3cab3cd1233bdcdf1785f4692c143a5be5fde7963198cfa2d74fc7d9c606", 
"2421e8fa7063ed95cfc14a295e2254e14aa5ca1e8fba117c5f4f27dad65cef18", 
"26f6e597bba2f2f1cc4e0d505e0fe7f16ad77271fd911b582b6350a4510a2eb1", 
"292bf9d2c3ee37eb4d413f23f12507e5871d47548cb72dff9f1a211b10679e32", 
"29b4c8366efb3e73bc76ffa8aa52013eb8160c98a18e1997aa906a05a15c4786", 
"2aeaebb41eaabceca08746649ab149f5c68b4b5b96a281159f10e0a967bb39af", 
"2ce4ba1cf5aa79249648246ffdfe46fdd8a42be8ec368dfc7989a89616ef3c26", 
"2e0a7bda343db6684115241dcf018dc6b1d16d563614f9e45c330567bc128fa1", 
"2e254b89d0d17b898d8f05bd84c2b982b346ba6241d41413686a79ab7d8da1fe", 
"3166a179de48204cfebd1cd70ad0d1f0d42769b93d2d72a0ad25a740aaa8da9f", 
"3492ba15e342e0665e53bbddc3738c7d2836ba13be3f4e4e41d9128778cce6fe", 
"368ca29d3ff7b234eb4b866b396367d5bad9ae8e71ee8737b351d057bdf2e632", 
"36b56876b006273d34b566c67508c488e49b429742c647bacc3ce55069af9d66", 
"3a091a12d80b1a9dd1d916af41fd1c74c168a27c4e7ef74a85a8b97e6b1c57f3", 
"3a25c7e1145e49c50a440bedd7ad4273b180e14b34089caced81791bcb757ab8", 
"3afc4fe7fb1fdc9f96b2087d12e899ac80746a8e5afa51e309c5ee2cae102cc6", 
"3be4f2bd0c5d99dca82fa9ef55e1a5447fa4fb006e8be90215e5bd2c2f838cf6", 
"3c8bcf4ce37c348517d1e1ac7d986e6993f46bfc1ef7cdc0aeca5bce164a5ed1", 
"3ca6ed74da619be17d401d14d312c569156c6550894cbbdc3d338741f8559132", 
"3d8f675f1890ceffd21688dcc6ef98cec729213c79759e4c80a33cc7ed54f6f6", 
"3da296dc5ba0010d3df89941771309126309d349f599b6caf8a9bb705820706f", 
"3ddcc03d81232fd0126f78c5cbcbc17e6898ec064a9944d7d7b28bee7ca0d878", 
"409a41987274a15c204602c5ae18e42b9a1c4b74fa434ee7a0d228f1c9794f65", 
"416e87167b7df7aac0b1d19b44aaae6f374f8d506ec54364b0c97b074458a0a8", 
"44c4aee40a12fd4fce4aa87febd8dbd689f06bb26e7700c492d9f728ef624ee1", 
"479e86260a31a04b61f1d284e0541fd8096b38c22fd18a19f070efe7adaa179d", 
"47d5637afbc09a7387de64c4a3050ef9d21d23928149c54955a7f08b101748e8", 
"49ef0e2104adc5dcfd59ec8cf539bc3d0922cdc769869209b91555e8fec7d349", 
"4a567ac55bd296932c89913914a01dc16d3ed3ca419680cf859829b4ec5ec662", 
"4dc72b4b0c1e89d0b20894f67716d8a71dad092188c7627965f5d5314126972c", 
"4ddc8ad0c1f28a2a878474c7187d89e5ae02ac4bea0c12bcabae9901ca18e461", 
"4fcdd8befe66235c0e8792eb6980e2880b7abd0123fedafa75eaa47903852796", 
"5048772db5bb399ef14140215c71b7d060ddf0c368261ad260ecb96d5def802b", 
"50cf272a7ff1fbdadafd85e67cdcf979c0d1e22a6f5d505831cbce39a189fa74", 
"51af67533d59abea5602b34ad1764b1bcd44b5e4f510a7b7ea45482b25417632", 
"51f2ccf8d5621a7cc4048311d393d406cae4ef26307a1397acd0993517a02187", 
"5344250ae4bca7ff8164e53d9fbebd674ed2af30c217b344a569bf65a45a5d99", 
"53ca44ea202c15fb2dfe6ab815871436ac61e22b9264e6d236ac2841a9ff1bcb", 
"54766fa4be2866c8f8397bac771aafbfc7106db476ede2c0ddb7d648c7f8856a", 
"5757d73e17118586df3e6d60cbf2315ff1dc688d2685e5f07845c4b56f44b6ba", 
"5799dab6dd3325039744fe0b81ac128a3d697cc8d14586d43fc20a8d6049c636", 
"59281d00f707101af80165dd97a302d9b38bb6fa8a3d1afdd8e31d7d9c08c1b6", 
"5a6e4d1a6de8d1e4c43745406fb9c6689c9357d56d47fc4e1ccc826d4ab9001e", 
"5a92ed01346a9de008b7943d29d788c8bf95d57a5368ebb60629b0436cd24b7d", 
"5aa79e26e86d5d074f65aee2feffae7050f61b2daaa5b01920678181d576a194", 
"5aebfd0f6d18cb695982fd71e8705202e42c911a0df79a49c47b8fa08176a342", 
"5debf54063b68f64901b1559dcc0dec09670ec76aaaa858131eabc2392c40976", 
"63418ab30d9d557d742edd44583364ac4c5e5fc15dc1863d07bb366e519d46e3", 
"63b6b54d0c1484867a0a72b641204bd10b945455621e573994135cf6377c3389", 
"63f9dabe4f2e0f664868de6fe7f9ce2443db7443f24021caf1b2cda06801b7cc", 
"652505448c9c38076c682ecee925063cd45bbe4ad6ec2c96f2ed4a99c7f006c3", 
"659c340280b4ae53640409046d6ba0e142dca124f2bd6c638c8324cf2d1b0838", 
"66765342c59d7f1a7aa729931b0f2af72cc0b71358ea7cdd3e8f39ed2294132b", 
"674726d192962840463fa4e47b530d5bee584c28c7a3331e0acda8c3ad0ddc1a", 
"67b8ffacc05bb4cd3cf2585db0ba0bdefc13bea93cbb2bace4aa14379890fb46", 
"6b035600c102b263e898e7c25d55c4b14cb821f45370220611c0788231160b4f", 
"6b32287b7a1bb2b079c9a5946a3f337e4caba4a78022903a9af5b66848792b29", 
"716cf8d45a14e56ed0d7c63b85d60ec6baced3e7c43e6a189a8b18a2be3e9f49", 
"73df60635de4da05addd6ea3c8a0cf3be1ca410668a63ec057693057da367c93", 
"751a72c5f71efe03adec5375c23a2b91ed3cdbacb7c8466320809f853abb59e6", 
"76c1f705d7125d55924c2ac4dae07abe911d0ffb4f5ca2e1d27a2f2186f3a071", 
"77183427f0d48323dc6aaf1105877b92e3bd3d4e736bc77abd8e87cff266e0c3", 
"7c83acf859263509d2dbbe0f07cb90266cd8d184ff88c83c972684af20354284", 
"80026dac5038bb76d5c75f500aeb807d439ca6a333146e6d7d117c4781a34089", 
"82e69c3dd8fef335b86b4f4fbe25a1ac27d48e5ad0d00dcc8d33f2282bb38b25", 
"8319a9336b81f76e5942fa232aeda15b14beb16a8fe62ce9dda66eca86cd6b03", 
"859fd8521edb09fdb970e21bc3cee040735ac3b94c2a51c825b8798f8236fd54", 
"863eb19ffe85217ac8efc5dd0c32b4c6c28f639f9bad2c84796823b7a5fb2637", 
"8701cc9de4d6ab0afcebb1767dd0d4a434906f6f66092ec5aa764aaa4486cc8a", 
"882f02585925ee6aba9966556552b591760f11f61201b74254bb8ce3fa1a91aa", 
"88c9c1ae45655c4505bbf37d9f17ba81674693e616227618ccd25ecced2fc397", 
"8a3d27b60ce1b1575b5f0de528d28b8a9c1ef7db799fa30881d1ff360f03e59a", 
"8b0c786a1703274bbcfe1a9f20f24b4411a49d2a9677682a0328b4fa81d02ebd", 
"8b36fc2862259034358f3b8432a70f031295015b6975ed83a20d7d1d8ac2f681", 
"8b8905a4f27dfe7c31585af78b3edc76941f450b914587dcde88830dcabeefa8", 
"8c45df956595ccb4e7f145d2cd34c82f816bec8b16b3e619d89ef5dcc7c70b93", 
"8ce465fe751f82c9d7908adae36e1ad97bb1f33b11373b3e8dd50ae0cffa8de6", 
"8d61df117050a2d45dc4a8d541ee5d4ffaf2e012584ea6c1fbe2ea94ef65d04f", 
"8db2e83fa9536bfe12586b88f3e344b9543096adc22bb8963263c4c30e8f1d7e", 
"8fa7b87093d30ea3298e257f8517d54460521ae583b3e4aa554708f010215614", 
"8ff4881771d6191d46d9ea46503e1f8515b4c437ea67adcbf12142258c432dd9", 
"9202884dded8b014d48101bd97aa6d80cec743c2437e5d4496158b2d13a19db8", 
"9794134dbbadc4c9802f77001b5357269d3df4c290a565f2d0b108c45abfce65", 
"97ac376fb9099cef6fabf948ca563ce5e70dbbff820897b53689a284910b0fb4", 
"981993c0fe52ed2c38380732392e8894452c0dae55d8e11cd458ecc7b8a45f4d", 
"982db16ab13704550b10d584334df5914c4fe5309a10211f3c25672f7200cb87", 
"987a14457dbc6a9a397ab46c82005b718c7dd8b10a6c3d7d8dc300e79c30732c", 
"9b078aa545e9a862f895d372b25b2fd7ea5550217b83d45e0ad91b7b969db31f", 
"a0eb91edd5daf7b2576f8eb3109a60ccb1f0280a5afcc19659bd6547efc5979d", 
"a258265be5f7f2ab78c3bbb9946f3713188bd6ef844bf96f5eac5429679036ca", 
"a294f117b9c905203baabd164efa30d3e6a2929446899a52c426ad35f67c7eb7", 
"a2b192de8292bd2f65cbaed42be122fc36272b0228ef3bce8c5f33ea8ed437cb", 
"a5a2cf6332da3bd1cc79c936cbbdd3de9b0d576f17ff6fbdaaafe896e9e88ca2", 
"a67a66b3e385611d077166b652bff025a8a4ecc5569e34bb6856d3dca2919fa8", 
"a878587878c3de4df64fba2209420097a7db38d9d41537695255463bc428fe4e", 
"a8c7b71fcff0f08f34160f3296c27d59fbfb808a67c2ec437b34e9d967668b5b", 
"a9c7c1b0c03849e605533aee147c0c2c5c456641ece164f9ac69cf9aacfeefa1", 
"ab3ae8dc4b3ae0be7df1dfdb703b72ea21a5f2ae0e6b947088b0f62e10567d70", 
"ac794c50715c2b2685f5ae8ba0f4cbc3d953fa19a5c86c7ec54a2451d7236c71", 
"af42d9d175547bba724f0708e8ff5919499c88e5527368d70135b2c3e828cba5", 
"b005ddbf6bfb6efe4cc1aaac8d72fbd3c68881e6a026c715cf1be70f6c70cf05", 
"b22ce16967b23f36e05ff814819418b6d3c8712ab1cc8c8d3377f72bd0d6af5b", 
"b3f9c8efcb7f9c0979235c66474113954848ada89a5f28442a6023554aaf3f12", 
"b4bba02cc28643c9c7405e9810deee792d08453215c78f9010a1ce79edff2e71", 
"b4fb1ac1cf75428a3a68b12ba48170ecae66e5b55658b8178b999c5b29d66c57", 
"b5dde7c32a5ce6432d6beecf3f2487ff891a69f89227679e04bffd7b435ef091", 
"b612a5ea4d3dcc69487906f396eed3046f433e90005d44c36b055ff11f602646", 
"b612bbafdb837896522c3d865091d0f15d3e5d2ad938d003361d18077e48e0a6", 
"b8c5ed0b7f7910ed05899dc929f7e7e1063f03f7782e09b7b39103b3f5f85acb", 
"bc266e810557e3e8023c2d583437921cd51a04c02f0032c69cb147e583661dd6", 
"beb8ab9ac10c8f534da7c9828037a8b7447b1c347a9d54737505f95ce22ddf21", 
"c0969d9d9fdd99fa58c4f347d7888738fd8f7a8313ba541557f450bf7449dce9", 
"c14bd9b31f1380e4d4e8e89c1d9000b5b52f67265f003b320dfb067228217b2d", 
"c3b5a8c1f1a72937047a741f53e0ea82054716c58d456b8a31d446bf757f672c", 
"c6c9188f2e68dfb65abad9b8e3a9715b4f067eb038e9350c6988b6ed49b49c66", 
"c72545e70eae8df1378e03264a71dc31ed17e693be7c1ae970865154ceeb6a15", 
"ca81226f53a7cdcc92fe4658a4ae3fafceebdfe17d7dcaf56fc1b992d2d8af53", 
"cabf2cbc8e1235c8d1885345cc9735263ab91531f2ca9a114a5e29fa7723bb5b", 
"cc2ccaba0be871525c840b58ef6cf8978ebeb113fada80af4e4f2bc24b8a9faf", 
"d05204d028979e7fce9420fb9157dd68152714ee1b54cf5ccaef0791461bc2b2", 
"d0f4b17e80dacf858a18f8751a560ca11076474e6da7f5ff1dcd023e8b9009ca", 
"d25db5c721797cb9b945dea5ef68fe8ca05c3f049801f1ba86706928a3afd992", 
"d639aafb7d0d19aed977cda8aaffa1544f7d9d246fd4855edf3cf4e5c60e411c", 
"d71b4d52b5bba45afaffc84266cd7c7880699f1b4c2c5b4a597985401699eb65", 
"d856a4924f9ddf174c6838b49c5e655b245681df8974fa8aed9797a5d8f7dbac", 
"d875b9c8a8ca2c0a901e32fc88c2cf2d5bbabddcb92a3579a9d1112176c06aac", 
"d977bea926a9a1a7a5dba6bc5b103b182606fea484e46ee7a449b8924e3b66a6", 
"d9d7451829ec38ca26b30373c4627f5d1d1588f363063b588b7f56523599e748", 
"da7afe33b380af5b490f77635c67a1000117575411c5f214b01a5cafd5520d98", 
"db701e9598a812dada990afbd489e8c9fa5ec71fd06232cf46dd97a75acf5bcf", 
"db851fb46fce284c2bc92a3710e767e955978b314d43ba123aa61260e430dfd5", 
"de210239948e447535998711e3fed8ec4a0ad63e2e785822a79e9bd67f2c46ed", 
"de83dc2c6de928fb609c97afe55ab09a1c434818476a28a8c44550254799eea9", 
"e1e3160e88450628aa76c4753259f327242df5a2e5809f9709ff70a70ff7b35f", 
"e574bbb5b38a22517a61f2605bb75c90fd7273b683890cf13b6612c2feb620ac", 
"e6302300219b7cc58f914c383206e240f8e07b166c57a92b47daa9fac516a4f2", 
"eb13b3e02db36ff1d73dcea51120834b1807e48fc64b322fbb50573bb7d736a5", 
"eb589c54ee9b571f8ad0b77d9fc8271e3c91ad3238356936d691c827716625e5", 
"ec55dd923483772d1b47417fb846ff4d6d695a5bc66c300cb748492a1b63956d", 
"ec824852c20916ca206cff722c3b9f07d3c34ff46d5e29e233735348cc9ed8a2", 
"ed0ae1a8c793c39134dca6625ec44839349e669e43f7492973b67f79e68184f1", 
"eee431466a3f084d7121812e4d5365c545122b214fb8c4b3dffa21fcc0660d2c", 
"f3bf6d46b1c7bfb92019e7e12528614bcd16d74e1eb6c2d9b052fea961041a83", 
"f3e7c69e13e91fb61c68008ba3ae3296056945748821f8bff7b3ac3da668bb72", 
"f42baf781aced50f5e3af8c93b0024dfc67c9b64fc69a7477bdcbb973de38aeb", 
"fb568c71e4a454b73217e6b314aea4ff9e8f146801307e886722383b760f1ee3", 
"fc4fea2a86852f5056910d4dca50dcfb3b68c084ecacc9180879b788e599bcc3", 
"fcc69a0efcecaf2536490389a424da8ab2a6d1f95639b287419d66810a9b29f1", 
"fdda0db48ea9ef00687c3c8b5a6c77f25b05bfe7049cea2dd823a877c298f83b", 
"fe00a8a871cb190708a7f4b405c527e2668ae27c666c0b57de975c286e252700"
)


ppoo_clus_pre_pam_om2_m_clasificaciones_sin_asw_neg<-
df_filled2[, c("run", "glosa_pueblo_originario")]|>
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2")], by="run", multiple="first")|> 
  dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "ci_cya_conadi", "ppoo_minsal_y_rsh_2010")], by="run", multiple="first")|>
  mutate(glosa_pueblo_originario_rec = if_else(
 glosa_pueblo_originario == "NINGUNO" & (ci_cya_conadi ==1 | ppoo_minsal_y_rsh_2010==1), "DESCONOCIDO", glosa_pueblo_originario
  )) |>
  filter(!run %in% run_aws_neg_m_c2) |> 
  group_by(run, clus_pam_om2) |>
  summarise(
    tiene_otro_po = any(grepl("AYMARA|COLLA|DIAGUITA|KAWÉSQAR|RAPA NUI|YAGÁN|OTRO", glosa_pueblo_originario_rec)),
    tiene_mapuche = any(grepl("MAPUCHE", glosa_pueblo_originario_rec)),
    tiene_ninguno = any(glosa_pueblo_originario_rec == "NINGUNO"),
    tiene_desconocido = any(glosa_pueblo_originario_rec == "DESCONOCIDO"),
    .groups = "drop"
  ) |>
  mutate(clas_por_run = case_when(
    tiene_mapuche ~ "Mapuche",
    !tiene_mapuche & tiene_otro_po ~ "Otros reportados (no mapuche)",
    tiene_ninguno & !tiene_otro_po & !tiene_mapuche ~ "Ninguno puro",
    tiene_desconocido & !tiene_ninguno & !tiene_otro_po & !tiene_mapuche ~ "Desconocido puro", T~NA_character_
  )) |>
  janitor::tabyl(clas_por_run, clus_pam_om2, show_na = T)

janitor::adorn_percentages(ppoo_clus_pre_pam_om2_m_clasificaciones_sin_asw_neg, "row") |> 
  mutate_if(is.numeric, ~round(.*100, 1)) |> 
  knitr::kable("markdown", caption="Otra forma de clasificación PPOO con reporte vs. conglomerados (por fila) (excluye ASW neg)")

cat("Fisher test (menos los -ASWs)\n")
janitor::fisher.test(ppoo_clus_pre_pam_om2_m_clasificaciones_sin_asw_neg, simulate.p.value=T, B=1e5)
#p-value = 0.00086

invisible("The Bonferroni corrected significant level for Exact approach is: 0.05/(R*C)=0.00625")
invisible("Para Otros reportados (no mapuche) p-value 0.01 en basado en residuo estandarizado, ajustado, y en la aproximación exacta, menor en TUS")
invisible("Para Mapuche p-value 0.03-0.05 en basado en aproximación exacta, residuo ajustado y estandarizado, respectivamente, de que es más alto proporción en TUS")
#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:

ppoo_clus_pre_pam_om2_m<-
  df_filled2[,c("run","glosa_pueblo_originario")]|> 
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om2","factor_inclusivo_real_hist_mas_autperc")], by="run", multiple="first")|> 
  dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario=="NINGUNO" & factor_inclusivo_real_hist_mas_autperc!="00"~ "DESCONOCIDO", T~glosa_pueblo_originario))|> 
  janitor::tabyl(glosa_pueblo_originario_rec, clus_pam_om2)|> 
  janitor::adorn_percentages("row")
#scale_fill_manual(values = rev(c("#D2B48C", "#E27A5B", "#708090", "#6B8E23", "#506070"  , "#2F4F4F", "#20B2AA"))) +
reshape2::melt(ppoo_clus_pre_pam_om2_m, id.vars = "glosa_pueblo_originario_rec")|> 
  dplyr::mutate(glosa_pueblo_originario_rec= 
    dplyr::recode(glosa_pueblo_originario_rec, 
      "OTRO (ESPECIFICAR)"="OTRO(n=77)", 
      "RAPA NUI (PASCUENSE)"="RAPA NUI(n=34)", 
      "YAGÁN (YÁMANA)"="YAGÁN(n=2)",
      "AYMARA"="AYMARA(n=13)",
      "COLLA"="COLLA(n=6)",
      "DIAGUITA"="DIAGUITA(n=3)",
      "KAWÉSQAR"="KAWÉSQAR(n=4)",
      "MAPUCHE"="MAPUCHE(n=255)",
      "DESCONOCIDO"=".DESCONOCIDO(n=1.985)",
      "NINGUNO"=".NINGUNO(n=9.156)"))|> 
ggplot(aes(x = glosa_pueblo_originario_rec, y = value, fill = variable)) + 
  geom_bar(stat = "identity", position = "fill") + 
  scale_fill_manual(values =  c(#, al menos un TSM(2)
  "6623, Primer mes, TSM y Comorbilidad(2)" = "#D2B48C", 
  "6612, Primer mes TUS(1)" = "#E27A5B"
)) +
  labs(title = NULL,
       x = "Grupo Étnico",
       y = "Proporción de Casos",
       fill = "Grupos") +  # Cambia el título de la leyenda a "Grupos"
  theme_minimal() +
  theme(
    axis.text.y = element_text(size = 12),           # Tamaño de las etiquetas de los grupos étnicos
    axis.text.x = element_text(size = 12),           # Tamaño de las etiquetas del eje X
    axis.title.x = element_text(size = 14),          # Tamaño del título del eje X
    axis.title.y = element_text(size = 14),          # Tamaño del título del eje Y
    plot.title = NULL,  # Tamaño y estilo del título del gráfico
    legend.title = element_text(size = 14, margin = margin(b = -.1)),          # Tamaño del título de la leyenda
    legend.spacing.y = unit(1.5, "lines"),
    legend.box.spacing = unit(0.5, "lines"),      # Controla el espacio entre la leyenda y el gráfico
    legend.margin = margin(5, 5, 5, 5),  
    legend.key.height = unit(1, "cm"),  
    legend.text = element_text(size = 12)            # Tamaño del texto de la leyenda
    ) + 
  coord_flip()  # Hacer el gráfico horizontal
ggsave("_figs/grafico_ancho_achatado_pam_om2_m_25.png", width = 10, height = 5, dpi=1000)

Sólo RAPA NUI y MAPUCHE (ver por RUN, al menos una observación que se identifique RAPA NUI o MAPUCHE)

    Pearson's Chi-squared test with Yates' continuity correction

data:  (ppoo_clus_pre_pam_om2_m_act)
X-squared = 4.2837, df = 1, p-value = 0.03848

Sólo RAPA NUI (ver por RUN, al menos una observación que se identifique RAPA NUI)

    Fisher's Exact Test for Count Data

data:  (ppoo_clus_pre_pam_om2_m_rapa)
p-value = 0.3708
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
 0.1777951 8.2063651
sample estimates:
odds ratio 
  1.694165 

son 11 casos!
Clasificando por Mapuche, otro no-mapuche, ninguno puro y desconocido puro

Porcentajes por columna, conglomerado vs. PPOO

clas_por_run	6623, Un trimestre, TSM(4)	6612, Un trimestre, TUS(3)	6522, Un semestre TSM(1)	6574, Comorbilidad un trimestre(2)
Ninguno puro	80.2	77.7	82.0	83.4
Desconocido puro	15.5	18.0	12.5	13.1
Otros reportados (no mapuche)	1.8	0.5	2.5	0.9
Mapuche	2.5	3.7	3.0	2.6

PPOO por cluster

    Fisher's Exact Test for Count Data with simulated p-value (based on
    1e+05 replicates)

data:  ppoo_clus_pre_pam_om2_m_clasificaciones
p-value = 0.00102
alternative hypothesis: two.sided

Comparación post-hoc, conglomerado-PPOO autorreporte

Dimension	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
Ninguno puro	1.76 (p=0.632)	-1.76 (p=0.632)
Desconocido puro	-1.97 (p=0.393)	1.97 (p=0.393)
Otros reportados (no mapuche)	2.68 (p=0.059)	-2.68 (p=0.059)
Mapuche	-2.08 (p=0.301)	2.08 (p=0.301)

PPOO por cluster

Tampoco resulta significativo en https://adaptivedesignstrial.shinyapps.io/posthoc/
Ahora descartando ASWs negativos

Otra forma de clasificación PPOO con reporte vs. conglomerados (por fila) (excluye ASW neg)

clas_por_run	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
Desconocido puro	86.3	13.7
Mapuche	82.8	17.2
Ninguno puro	88.4	11.6
Otros reportados (no mapuche)	96.4	3.6

PPOO por cluster

Fisher test (menos los -ASWs)

    Fisher's Exact Test for Count Data with simulated p-value (based on
    1e+05 replicates)

data:  ppoo_clus_pre_pam_om2_m_clasificaciones_sin_asw_neg
p-value = 8e-04
alternative hypothesis: two.sided

PPOO por cluster

2.1.1. Trayectorias

Vemos los gráficos de las trayectorias

categories_pam_om2_m<-attr(States_Wide.seq_month_t_prim_adm_cens, "labels")
new_labels <- categories_pam_om2_m
new_labels[which(categories_pam_om2_m == "Otras causas")] <- "Otras\ncausas"
#new_labels[which(categories == "Consumo\nde sustancias")] <- "Consumo de\nsustancias"

# Creamos un vector con las columnas llenando con NA si faltan valores
sil_pam_om_clus2_m <- wcSilhouetteObs(as.dist(dist_month_om), 
        pamRange_month_om$clustering$cluster2, measure="ASW")


seq_plot_pam_om2_m <- ggseqplot::ggseqiplot(States_Wide.seq_month_t_prim_adm_cens, 
                                 group= ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2,
                                 facet_ncol=1, facet_nrow=2, sortv=sil_pam_om_clus2_m) +
  theme(legend.position = "none")+
  labs(x="Meses", y="# IDs de usuarios")+
  #guides(fill = guide_legend(nrow = 1))+
  theme(
    panel.spacing = unit(0.1, "lines"),  # Reduce el espaciado entre los paneles
    axis.text.y = element_text(size = 15),           # Tamaño de las etiquetas de los grupos étnicos
    axis.text.x = element_text(size = 15),           # Tamaño de las etiquetas del eje X
    axis.title.x = element_text(size = 15),          # Tamaño del título del eje X
    axis.title.y = element_text(size = 15, margin = margin(r = -10)),#,margin = margin(l = -10)),
    strip.text = element_text(size = 15, margin = margin(b =-15, t=10)),
    legend.text = element_text(size = 15),
    legend.spacing.x = unit(0.1, 'cm'),  # Alinea el título de la leyenda hacia la izquierda
    legend.box.margin = margin(t = 0, r = 0, b = 0, l = -50),
    legend.position = "bottom", 
    legend.justification = "left",
    panel.spacing.y = unit(0.2, "lines"),
    plot.margin = margin(10, 10, 10, 10), # Ajusta márgenes globales
    strip.placement = "outside",   # Para colocar las tiras fuera de los ejes
    strip.background = element_blank() # Elimina el fondo para que parezca más espacioso
    #legend.key.size = unit(1.5, "lines"),        # Aumenta el tamaño de los símbolos en la leyenda
  )+
  guides(fill = guide_legend(nrow = 1)) +
  scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513",
                                                  "#FFFFFF","#808080","#000000"))+
  scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513",
                                                   "#FFFFFF","#808080","#000000"))
seq_plot_pam_om2_m 

ggsave(filename="_figs/clusters_pam_om2_m_mod_25.png", seq_plot_pam_om2_m,  width = 11, height = 5.5, dpi=1000)

Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)

seq_plot2_pam_om2_m <- ggseqdplot(States_Wide.seq_month_t_prim_adm_cens, 
                                 group= ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2,
                                 facet_ncol=1, facet_nrow=2) +
  theme(legend.position = "none")+  # Colocar la leyenda abajo
  labs(x="Meses", y="Frecuencia relativa de estados")+
  theme(
    panel.spacing = unit(0.1, "lines"),
    axis.text.y = element_text(size = 15),           # Tamaño de las etiquetas de los grupos étnicos
    axis.text.x = element_text(size = 15),           # Tamaño de las etiquetas del eje X
    axis.title.x = element_text(size = 15),          # Tamaño del título del eje X
    axis.title.y = element_text(size = 15, margin = margin(r = -5)),
    strip.text = element_text(size = 15),
    panel.spacing.y = unit(0.5, "lines"),
    strip.placement = "outside",   # Para colocar las tiras fuera de los ejes
    strip.background = element_blank() # Elimina el fondo para que parezca más espacioso
    #legend.key.size = unit(1.5, "lines"),        # Aumenta el tamaño de los símbolos en la leyenda    
  )  # Colocar la leyenda abajo
seq_plot2_pam_om2_m
ggsave("_figs/clusterspam_om22_m_mod_25.png",seq_plot2_pam_om2_m, width = 11, height = 5.5, dpi=1000)


table_data_pam_om2_m <- sprintf("%1.2f",pamRange_month_om$stats[1,])
table_data_pam_om2_m <-as.data.frame(t(table_data_pam_om2_m))
colnames(table_data_pam_om2_m)<-attr(pamRange_month_om$stats, "name")
table_data_pam_om2_m|> knitr::kable()

PBC	HG	HGSD	ASW	ASWw	CH	R2	CHsq	R2sq	HC
0.16	0.38	0.38	0.40	0.40	833.05	0.11	769.36	0.10	0.49

Trayectorias de hospitalización, frecuencia relativa de estados en un gráfico de barras apiladas por trimestre.

Trayectorias de hospitalización, frecuencia relativa de estados en un gráfico de barras apiladas por trimestre.

De este modo, presenta el cambio agregado en la distribución de estados a lo largo del tiempo, sin considerar las secuencias individuales.

invisible("Definimos las observaciones que tienen siluetas negativas")
sil_neg_pam_om_clus2_m <- which(sil_pam_om_clus2_m<0)

invisible("A qué conglomerados pertenecen?")
table(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[sil_neg_pam_om_clus2_m, "clus_pam_om2"])

clus_pam_om2
6623, Primer mes, TSM y Comorbilidad(2)                 6612, Primer mes TUS(1) 
                                    171                                       0 

2.1.2.Exploración transiciones

2.1.2.a Transiciones- RM y no RM

Tasas de transición no RM a RM y viceversa

invisible("Tasas de transición no RM a RM y viceversa")

trim_tasa_pam_om2_m_cens_cnt<-  
  seqcount_t(States_Wide.seq_month_t_prim_adm_RM_cens, 
             group=ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2)|> 
  dplyr::filter(count>0)|> 
  dplyr::mutate(trans = paste0(from,"_", to))|> 
  dplyr::mutate(across(c("from","to"),~  gsub("\\[->\\s*|\\s*->\\s*\\]|\\[|\\]", "", .))) 
trim_tasa_pam_om2_m_cens_rate<-  
  seqtrate_t(States_Wide.seq_month_t_prim_adm_RM_cens, 
             group=ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2)|> 
  dplyr::filter(rate>0)|> 
  dplyr::mutate(trans = paste0(from,"_", to))|> 
  dplyr::mutate(across(c("from","to"),~  gsub("\\[->\\s*|\\s*->\\s*\\]|\\[|\\]", "", .)))

trim_tasa_pam_om2_m_cens_rate|>   
  dplyr::left_join(trim_tasa_pam_om2_m_cens_cnt, by=c("from"="from", "glosa_sexo"="glosa_sexo","to"="to"))|> 
  dplyr::rename("recuento"="count")|> 
  dplyr::filter(from %in% c("RM", "noRM"))|>  
  ggplot(aes(x = from, y = to, fill = rate, size=log(recuento+1))) +
  geom_tile() +
  coord_flip()+
  scale_fill_gradient(low = "white", high = "blue") + # Ajusta la escala de colores según tus preferencias
  labs(title = "Tasas de transición, Mes (s/censura)",
       x = "Desde",
       y = "Hacia",
       fill = "Rate") +
  theme_minimal() +
  facet_wrap(~glosa_sexo)+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))+
  geom_text(aes(label = sprintf("%1.2f", rate), size =log(recuento+1)*.5), color = "black")

invisible("Hay muy pocos casos que se entrecruzan entre noRM y RM (fuera de la diagnonal)")

Porcentajes de transición no-RM y RM por cada cluster

Hay muy pocos casos que se entrecruzan entre noRM y RM (fuera de la diagnonal)

2.1.2.b Transiciones

trim_tasa2_pam_om2_m_cens_cnt<-  
   seqcount_t(States_Wide.seq_month_t_prim_adm_cens, 
             group=ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2)|> 
  dplyr::filter(count>0)|> 
  dplyr::mutate(trans = paste0(from,"_", to))|> 
  dplyr::mutate(across(c("from","to"),~  gsub("\\[->\\s*|\\s*->\\s*\\]|\\[|\\]", "", .))) 
trim_tasa2_pam_om2_m_cens_rate<-  
  seqtrate_t(States_Wide.seq_month_t_prim_adm_cens, 
             group=ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2)|> 
  dplyr::filter(rate>0)|> 
  dplyr::mutate(trans = paste0(from,"_", to))|> 
  dplyr::mutate(across(c("from","to"),~  gsub("\\[->\\s*|\\s*->\\s*\\]|\\[|\\]", "", .)))

trim_tasa2_pam_om2_m_cens_rate|>   
  dplyr::left_join(trim_tasa2_pam_om2_m_cens_cnt, by=c("from"="from", "glosa_sexo"="glosa_sexo","to"="to"))|> 
  dplyr::rename("recuento"="count")|> 
  #dplyr::filter(from %in% c("RM", "noRM"))|>  
  ggplot(aes(x = from, y = to, fill = rate, size=log(recuento+1))) +
  geom_tile() +
  coord_flip()+
  scale_fill_gradient(low = "white", high = "blue") + # Ajusta la escala de colores según tus preferencias
  labs(title = "Tasas de transición, Mes (s/censura)",
       x = "Desde",
       y = "Hacia",
       fill = "Rate") +
  theme_minimal() +
  facet_wrap(~glosa_sexo)+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))+
  geom_text(aes(label = sprintf("%1.2f", rate), size =log(recuento+1)*.5), color = "black")

Porcentajes de transición, transiciones posteriores, por cada cluster

2.1.2.c Tiempo promedio por cluster

seq_mean_t(States_Wide.seq_month_t_prim_adm_cens, 
             group=ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2)|> 
  data.table::as.data.table(keep.rowname=T)|>
  dplyr::mutate(rn= gsub("\\d", "", rn)) |> 
  ggplot(aes(x=rn, fill= factor_inclusivo, y=Mean))+
  geom_bar(width = 1, stat = "identity") +
  theme_minimal() +
  facet_wrap(~factor_inclusivo)+
  labs(title = NULL,
       x = NULL,
       y = NULL) +
  scale_fill_manual(values = rev(c("#D2B48C", "#E27A5B"))) +
  coord_flip()+
  theme(#axis.text.x = element_blank(),
    #axis.text.y = element_blank(),
    panel.grid = element_blank()) +
#  scale_fill_brewer(palette = "Pastel1", labels=c("Sin\nautoidentificación\nni reconocimiento", "Autoidentificación\nsin reconocimiento", "Ambas")) +
  geom_text(aes(label = round(Mean,1)), 
            position = position_stack(vjust = 0.5), 
            size = 3.5, # Ajusta el tamaño de la fuente aquí
            color = "black", # Color del texto
            family = "sans", # Puedes cambiar la fuente si lo deseas
            background = element_rect(fill = "white", color = NA)) + # Fondo blanco
  theme(legend.title = element_blank(), legend.position="none")

Tiempo promedio en cada estado por estatus PPOO (Trimestral c/censura)

Observamos que aquellos en el conglomerado que se encuentra ingresada por trastornos de salud mental tiene un porcentaje levemente mayor de censura.

2.1.3. Comparación variables

2.1.3.a. Comparación covariables- PPOO

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|>
  dplyr::count(clus_pam_om2, factor_inclusivo_real_hist_mas_autperc)|>
  dplyr::group_by(clus_pam_om2)|>
  dplyr::mutate(n_prop = paste0(n, " (",scales::percent(n / sum(n), accuracy=.1),")"))|>
  dplyr::select(-n)|>
  tidyr::pivot_wider(names_from = factor_inclusivo_real_hist_mas_autperc, values_from = n_prop, values_fill = "0")|> 
  knitr::kable("markdown", col.names=c("Conglomerados","No se identifica/no pertenece", "No se identifica/hay reconocimiento", "Se identifica/hay reconocimeinto"), caption="Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO") # leve mayor identificación y reconocimiento en TUS

Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO

Conglomerados	No se identifica/no pertenece	No se identifica/hay reconocimiento	Se identifica/hay reconocimeinto
6623, Primer mes, TSM y Comorbilidad(2)	4712 (80.5%)	661 (11.3%)	484 (8.3%)
6612, Primer mes TUS(1)	598 (77.8%)	93 (12.1%)	78 (10.1%)

Vemos las categorías de clasificación de PPOO según autopercepción (en MINSAL y en RSH) y reconocimiento CONADI.

chisq_cramerv(table(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$factor_inclusivo_real_hist_mas_autperc,ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2))
# $chisq_statistic
# [1] "3.83"
# 
# $chisq_df
# df 
#  2 
# 
# $chisq_p_value
# [1] "0.1471"
# 
# $cramers_v
# [1] "0.02"

cat("Descartando valores negativos en sil width")
chisq_cramerv(
with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m) , table(factor_inclusivo_real_hist_mas_autperc, clus_pam_om2))
)
#[1] "0.1532"

$chisq_statistic
[1] "3.83"

$chisq_df
df 
 2 

$chisq_p_value
[1] "0.1471"

$cramers_v
[1] "0.02"

Descartando valores negativos en sil width$chisq_statistic
[1] "3.75"

$chisq_df
df 
 2 

$chisq_p_value
[1] "0.1532"

$cramers_v
[1] "0.02"

No se constata una asocaición entre autoidentificación/reconocimiento perteneciente a PPOO y la pertenencia a conglomerados.

Ahora lo hacemos con la versión binarizada

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$factor_inclusivo_real_hist_mas_autperc_bin<- ifelse(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$factor_inclusivo_real_hist_mas_autperc=="00",0,1)

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|>
    dplyr::count(clus_pam_om2, factor_inclusivo_real_hist_mas_autperc_bin)|>
    dplyr::group_by(clus_pam_om2)|>
    dplyr::mutate(n_prop = paste0(n, " (",scales::percent(n / sum(n), accuracy=0.1),")"))|>
    dplyr::select(-n)|>
    tidyr::pivot_wider(names_from = factor_inclusivo_real_hist_mas_autperc_bin, values_from = n_prop, values_fill = "0")|>  #rio::export("clipboard")
  knitr::kable("markdown", caption="Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO / Binario")

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
    dplyr::mutate(death_time_rec=ifelse(death_time==60,0,1))|> 
    janitor::tabyl(factor_inclusivo_real_hist_mas_autperc_bin,clus_pam_om2)|>
    janitor::chisq.test()
#X-squared = 2.9179, df = 1, p-value = 0.0876
#
#
cat("Descartando valores negativos en sil width")
chisq_cramerv(
  with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m) , table(factor_inclusivo_real_hist_mas_autperc_bin, clus_pam_om2))
)# [1] "0.0938"

Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO / Binario

clus_pam_om2	0	1
6623, Primer mes, TSM y Comorbilidad(2)	4712 (80.5%)	1145 (19.5%)
6612, Primer mes TUS(1)	598 (77.8%)	171 (22.2%)

    Pearson's Chi-squared test with Yates' continuity correction

data:  janitor::tabyl(dplyr::mutate(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens,     death_time_rec = ifelse(death_time == 60, 0, 1)), factor_inclusivo_real_hist_mas_autperc_bin,     clus_pam_om2)
X-squared = 2.9179, df = 1, p-value = 0.0876

Descartando valores negativos en sil width$chisq_statistic
[1] "2.81"

$chisq_df
df 
 1 

$chisq_p_value
[1] "0.0938"

$cramers_v
[1] "0.02"

Tampoco se observa asociación alguna. Hicimos una prueba post-hoc usando Bonferroni

2.1.3.b. Comparación covariables- Mortalidad

# 
invisible("No hay nada, el tiempo promedio de censura es similar")

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
  dplyr::mutate(death_time_rec=ifelse(death_time==60,0,1))|> 
  janitor::tabyl(clus_pam_om2,death_time_rec)|> 
  dplyr::mutate(`1`= paste0(`1`," (", scales::percent(`1`/(`0`+`1`), accuracy=.1),")"))|> 
  dplyr::left_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
  dplyr::group_by(clus_pam_om2)|> 
  dplyr::summarise(mean=sprintf("%1.1f",mean(cens_time)),
  median=sprintf("%1.1f",quantile(cens_time, .5)), p25=sprintf("%1.1f",quantile(cens_time, .25)), p75=sprintf("%1.1f",quantile(cens_time, .75))), by="clus_pam_om2")|> 
  dplyr::select(-`0`)|> 
  knitr::kable("markdown", col.names=c("Conglomerado","Mortalidad observada", "Promedio", "Mediana", "Q1", "Q3"), caption="Post-hoc, conglomerado vs. Mortalidad y tiempo a censura")

Post-hoc, conglomerado vs. Mortalidad y tiempo a censura

Conglomerado	Mortalidad observada	Promedio	Mediana	Q1	Q3
6623, Primer mes, TSM y Comorbilidad(2)	78 (1.3%)	53.8	53.8	50.9	56.5
6612, Primer mes TUS(1)	21 (2.7%)	54.6	54.9	51.6	57.8

# Cargar las librerías necesarias
library(survival)
library(ggplot2)

# Crear la variable de supervivencia
surv_obj <- Surv(time = ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$death_time,
                 event = ifelse(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$death_time==60,0,1))

cat("sin siluetas negativas")

sin siluetas negativas

# Realizar el análisis de Log-Rank (survdiff)
survdiff(surv_obj ~ clus_pam_om2,
                      data = ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens)

Call:
survdiff(formula = surv_obj ~ clus_pam_om2, data = ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens)

                                                        N Observed Expected
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2) 5857       78     87.6
clus_pam_om2=6612, Primer mes TUS(1)                  769       21     11.4
                                                     (O-E)^2/E (O-E)^2/V
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2)      1.05      9.09
clus_pam_om2=6612, Primer mes TUS(1)                      8.04      9.09

 Chisq= 9.1  on 1 degrees of freedom, p= 0.003 

# Chisq= 9.1  on 1 degrees of freedom, p= 0.003 
 survdiff(Surv(death_time, ifelse(death_time==60,0,1)) ~ clus_pam_om2,
                      data = slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m) )

Call:
survdiff(formula = Surv(death_time, ifelse(death_time == 60, 
    0, 1)) ~ clus_pam_om2, data = slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, 
    -sil_neg_pam_om_clus2_m))

                                                        N Observed Expected
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2) 5686       76     85.5
clus_pam_om2=6612, Primer mes TUS(1)                  769       21     11.5
                                                     (O-E)^2/E (O-E)^2/V
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2)      1.06      8.95
clus_pam_om2=6612, Primer mes TUS(1)                      7.89      8.95

 Chisq= 8.9  on 1 degrees of freedom, p= 0.003 

# Chisq= 8.9  on 1 degrees of freedom, p= 0.003 
 
# Ajustar el modelo de Kaplan-Meier
km_fit <- survfit(surv_obj ~ clus_pam_om2,
                  data = ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens)

# Extraer los datos del modelo Kaplan-Meier para usar con ggplot
km_data_m <- data.frame(
  time = km_fit$time,
  surv = km_fit$surv,
  upper = km_fit$upper,
  lower = km_fit$lower,
  strata = rep(c("6623, Primer mes, TSM y Comorbilidad(2)","6612, Primer mes TUS(1)"), km_fit$strata)
)

biostat3::survRate(Surv(time = death_time,
                        event = ifelse(death_time==60,0,1)) ~ clus_pam_om2, 
                   data= ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens) |> 
    dplyr::mutate(across(c("rate", "lower", "upper"),~sprintf("%1.2f",.*10000))) 

                                                                                clus_pam_om2
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2) 6623, Primer mes, TSM y Comorbilidad(2)
clus_pam_om2=6612, Primer mes TUS(1)                                 6612, Primer mes TUS(1)
                                                         tstop event rate lower
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2) 348073.88    78 2.24  1.77
clus_pam_om2=6612, Primer mes TUS(1)                  45240.02    21 4.64  2.87
                                                     upper
clus_pam_om2=6623, Primer mes, TSM y Comorbilidad(2)  2.80
clus_pam_om2=6612, Primer mes TUS(1)                  7.10

epitools::rateratio.midp(c(21 , 78,  45240.02 , 348073.88))

$data
          Outcome
Predictor  Cases Person-time
  Exposed1    21    45240.02
  Exposed2    78   348073.88
  Total       99   393313.90

$measure
          rate ratio with 95% C.I.
Predictor   estimate     lower     upper
  Exposed1 1.0000000        NA        NA
  Exposed2 0.4799882 0.3020483 0.7982259

$p.value
          two-sided
Predictor   midp.exact       wald
  Exposed1          NA         NA
  Exposed2 0.005725614 0.00246084

attr(,"method")
[1] "Median unbiased estimate & mid-p exact CI"

# Crear el gráfico de Kaplan-Meier con ggplot2
ggplot(km_data_m, aes(x = time, y = surv, color = strata)) +
  geom_step(size = 1.2) +  # Curvas de supervivencia
  geom_ribbon(aes(ymin = lower, ymax = upper, fill = strata), alpha = 0.2, color = NA) +  # Intervalos de confianza
  labs(
    title = "Curvas de Kaplan-Meier",
    x = "Tiempo (meses)",
    y = "Probabilidad de Supervivencia",
    color = "Grupo",
    fill = "Grupo"
  ) +
  theme_minimal() +
  theme(legend.position = "bottom") +
  scale_color_manual(values = c("#E2725B", "#D2B48C")) +  # Colores para las curvas
  scale_fill_manual(values = c("#E2725B", "#D2B48C"))    # Colores para las áreas sombreadas

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
             dplyr::mutate(death_time_rec=ifelse(death_time==60,0,1))|> 
             janitor::tabyl(death_time_rec,clus_pam_om2)|> 
              janitor::chisq.test(correct=T)
#X-squared = 8.1148, df = 1, p-value = 0.00439

ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
             dplyr::mutate(death_time_rec=ifelse(death_time==60,0,1))|> 
             janitor::tabyl(death_time_rec,clus_pam_om2)|> 
              janitor::fisher.test(simulate.p.value=T, B=1e5)
# odds ratio 
#   2.079773 
#p-value = 0.006252

cat("Descartando valores negativos en sil width")
chisq_cramerv(
with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m)|> 
             dplyr::mutate(death_time_rec=ifelse(death_time==60,0,1)), table(death_time_rec , clus_pam_om2))
)
# $chisq_statistic
# [1] "7.98"
# 
# $chisq_df
# df 
#  1 
# 
# $chisq_p_value
# [1] "0.0047"
# 
# $cramers_v
# [1] "0.04"

#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
##_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
tab_cl_mortalidad_pam_om2_m<- 
  ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
  dplyr::mutate(death_time_rec=ifelse(death_time==60,0,1)) |> 
  janitor::tabyl(death_time_rec,clus_pam_om2) |> 
  as.matrix(ncol=2)

labels_pam_om2_m <- c(
  "6623, Primer mes, TSM y Comorbilidad(2)",
  "6612, Primer mes TUS(1)"
)

    Pearson's Chi-squared test with Yates' continuity correction

data:  janitor::tabyl(dplyr::mutate(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens,     death_time_rec = ifelse(death_time == 60, 0, 1)), death_time_rec,     clus_pam_om2)
X-squared = 8.1148, df = 1, p-value = 0.00439


    Fisher's Exact Test for Count Data

data:  janitor::tabyl(dplyr::mutate(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens,     death_time_rec = ifelse(death_time == 60, 0, 1)), death_time_rec,     clus_pam_om2)
p-value = 0.006252
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
 1.212205 3.426554
sample estimates:
odds ratio 
  2.079773 

Descartando valores negativos en sil width$chisq_statistic
[1] "7.98"

$chisq_df
df 
 1 

$chisq_p_value
[1] "0.0047"

$cramers_v
[1] "0.04"

Se constata una asociación muy modesta entre la pertenencia a un conglomerado y mortalidad, siendo 6025, Un trimestre, TUS(1) mayores sus chances de mortalidad.

2.1.3.c. Comparación covariables- no RM vs. RM

chisq_cramerv(with(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, 
  table(clus_pam_om2,codigo_region_rec_base)
  ))
# $chisq_statistic
# [1] "32.05"
# 
# $chisq_df
# df 
#  1 
# 
# $chisq_p_value
# [1] "<0.001"
# 
# $cramers_v
# [1] "0.07"

cat("Descartando valores negativos en sil width")
chisq_cramerv(
with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m), table(codigo_region_rec_base , clus_pam_om2))
)
# $chisq_p_value
# [1] "<0.001"

table(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2,ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$codigo_region_rec_base)|> 
    data.frame() |>  
    dplyr::group_by(Var1) |> 
    dplyr::mutate(perc= scales::percent(Freq/sum(Freq), accuracy=.1)) |> 
    dplyr::ungroup() |> 
    dplyr::mutate(Freq= Freq)  |>
    dplyr::mutate(fp= paste0(Freq, " (", perc,")" )) |> 
    dplyr::select(-Freq, -perc) |> 
    tidyr::pivot_wider(names_from="Var2",values_from="fp")  |> 
  knitr::kable("markdown", caption="Frecuencias absolutas y relativas por fila", col.names= c("Conglomerados", "Fuera de RM", "En RM"))

$chisq_statistic
[1] "32.05"

$chisq_df
df 
 1 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.07"

Descartando valores negativos en sil width$chisq_statistic
[1] "29.68"

$chisq_df
df 
 1 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.07"

Frecuencias absolutas y relativas por fila

Conglomerados	Fuera de RM	En RM
6623, Primer mes, TSM y Comorbilidad(2)	3063 (52.3%)	2794 (47.7%)
6612, Primer mes TUS(1)	486 (63.2%)	283 (36.8%)

2.1.3.e. Comparación covariables- Región

tab_cluster_region_pam_om2_m<-
ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
  dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
  dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], 
                    by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first")|> 
  janitor::tabyl(codigo_region, clus_pam_om2)|> 
  janitor::adorn_percentages("col")|> 
  janitor::adorn_rounding(digits = 2)

#colnames(tab_cluster_region_pam_om4_q)<- c("reg", "c1", "c4", "c3", "c5", "c6", "c7", "c8", "c9", "c2")
cod_reg_homo_pam_om2_m<-
data.frame(
  codigo_region = 1:16,
  nombre_region = c(
    "Región de Tarapacá",
    "Región de Antofagasta",
    "Región de Atacama",
    "Región de Coquimbo",
    "Región de Valparaíso",
    "Región del Libertador General Bernardo O'Higgins",
    "Región del Maule",
    "Región del Biobío",
    "Región de La Araucanía",
    "Región de Los Lagos",
    "Región de Aysén del General Carlos Ibáñez del Campo",
    "Región de Magallanes y de la Antártica Chilena",
    "Región Metropolitana de Santiago",
    "Región de Los Ríos",
    "Región de Arica y Parinacota",
    "Región de Ñuble"
  ),
  stringsAsFactors = FALSE
)

dplyr::mutate(tab_cluster_region_pam_om2_m, promedio_fila = rowMeans(across(2:length(colnames(tab_cluster_region_pam_om2_m)))))|> 
  dplyr::arrange(desc(promedio_fila))|> 
  dplyr::left_join(cod_reg_homo_pam_om2_m, by="codigo_region")|> 
  dplyr::select(codigo_region, nombre_region, everything())|> 
  dplyr::select(-promedio_fila)|> 
  dplyr::mutate_at(3:(length(colnames(tab_cluster_region_pam_om2_m))+1),~scales::percent(.))|> 
  knitr::kable(caption="Porcentaje por región")

Porcentaje por región

codigo_region	nombre_region	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
13	Región Metropolitana de Santiago	46%	35%
10	Región de Los Lagos	6%	19%
5	Región de Valparaíso	8%	13%
8	Región del Biobío	10%	8%
9	Región de La Araucanía	5%	4%
6	Región del Libertador General Bernardo O’Higgins	4%	3%
7	Región del Maule	4%	3%
14	Región de Los Ríos	3%	3%
1	Región de Tarapacá	2%	2%
11	Región de Aysén del General Carlos Ibáñez del Campo	2%	2%
12	Región de Magallanes y de la Antártica Chilena	2%	2%
16	Región de Ñuble	2%	2%
2	Región de Antofagasta	2%	1%
3	Región de Atacama	2%	1%
4	Región de Coquimbo	1%	1%
15	Región de Arica y Parinacota	2%	0%

tab_clus_reg_pam_om2_m<-
ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
    dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
    dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], 
                      by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first")|> 
    janitor::tabyl(codigo_region, clus_pam_om2)

chisq_cramerv(tab_clus_reg_pam_om2_m[,-1])
# $chisq_statistic
# [1] "184.78"
# 
# $chisq_df
# df 
# 15 
# 
# $chisq_p_value
# [1] "<0.001"
# 
# $cramers_v
# [1] "0.17"

janitor::fisher.test(tab_clus_reg_pam_om2_m, simulate.p.value=T, B=1e5)
#p-value = 1e-05

cat("Descartando valores negativos en sil width")
chisq_cramerv(
with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m)|> 
       dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
    dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first") , table(codigo_region , clus_pam_om2))
)#[1] "<0.001"

$chisq_statistic
[1] "208.18"

$chisq_df
df 
15 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.18"


    Fisher's Exact Test for Count Data with simulated p-value (based on
    1e+05 replicates)

data:  tab_clus_reg_pam_om2_m
p-value = 1e-05
alternative hypothesis: two.sided

Descartando valores negativos en sil width$chisq_statistic
[1] "210.80"

$chisq_df
df 
15 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.18"

Se encuentra una asociación significativa aunque débil entre región y conglomerados.

Por macrozona

tab_clus_macrozona_pam_om2_m<-
ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
  dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
  dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], 
                    by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first")|> 
dplyr::mutate(macrozona = dplyr::case_when(
  codigo_region %in% c(15, 1, 2, 3) ~ "Macrozona Norte",
  codigo_region %in% c(4, 5)~ "Macrozona Centro",
  codigo_region %in% c(6, 7, 16, 8) ~ "Macrozona Centro Sur",
  codigo_region %in% c(9, 14, 10) ~ "Macrozona Sur",
  codigo_region %in% c(11, 12) ~ "Macrozona Austral",
  TRUE ~ "RM"  # En caso de que algún código no esté especificado
))|> 
  janitor::tabyl(macrozona, clus_pam_om2) 

chisq_cramerv(dplyr::filter(tab_clus_macrozona_pam_om2_m,macrozona!="RM")[,-1])
# $chisq_statistic
# [1] "100.16"
# 
# $chisq_df
# df 
#  5 
# 
# $chisq_p_value
# [1] "<0.001"
# 
# $cramers_v
# [1] "0.13"

cat("Descartando valores negativos en sil width")
chisq_cramerv(
with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m)|> 
       dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
  dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], 
                    by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first")|> 
dplyr::mutate(macrozona = dplyr::case_when(
  codigo_region %in% c(15, 1, 2, 3) ~ "Macrozona Norte",
  codigo_region %in% c(4, 5)~ "Macrozona Centro",
  codigo_region %in% c(6, 7, 16, 8) ~ "Macrozona Centro Sur",
  codigo_region %in% c(9, 14, 10) ~ "Macrozona Sur",
  codigo_region %in% c(11, 12) ~ "Macrozona Austral",
  TRUE ~ "RM"  # En caso de que algún código no esté especificado
)) |> dplyr::filter(macrozona!="RM") , table(macrozona, clus_pam_om2))
)#[1] "<0.001"


cat("Descriptivos")
dplyr::filter(tab_clus_macrozona_pam_om2_m,macrozona!="RM")|>
  janitor::adorn_percentages("col")|>
  janitor::adorn_rounding(digits = 3)|>
  dplyr::mutate(across(-macrozona, ~ sprintf("%d (%.1f%%)", dplyr::filter(tab_clus_macrozona_pam_om2_m,macrozona!="RM")[[cur_column()]], . * 100))) |> 
  knitr::kable(caption="Porcentajes por columna, conglomerado vs. macrozona")

$chisq_statistic
[1] "78.50"

$chisq_df
df 
 4 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.15"

Descartando valores negativos en sil width$chisq_statistic
[1] "80.51"

$chisq_df
df 
 4 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.15"

Descriptivos

Porcentajes por columna, conglomerado vs. macrozona

macrozona	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
Macrozona Austral	197 (6.2%)	32 (6.4%)
Macrozona Centro	540 (17.0%)	105 (21.0%)
Macrozona Centro Sur	1173 (36.9%)	125 (25.1%)
Macrozona Norte	469 (14.7%)	34 (6.8%)
Macrozona Sur	804 (25.3%)	203 (40.7%)

chisq.posthoc.test(dplyr::filter(tab_clus_macrozona_pam_om2_m, macrozona != "RM")[-1])|>
  (\(df) {dplyr::mutate_at(df, 3:ncol(df), ~ round(as.numeric(gsub("\\*", "", .)), 3))})()|>
  group_by(Dimension)|>
  (\(df) {summarise(df, across(2:(ncol(df) - 1), ~paste0(first(sprintf("%1.2f", .)), " (p=", last(sprintf("%1.3f", .)), ")")))})()|>
  (\(df) {dplyr::mutate_at(df, 2:length(names(df)), ~gsub("p\\=0.000)", "p<0.001)", .))})()|>
  dplyr::mutate(Dimension= dplyr::filter(tab_clus_macrozona_pam_om6_q, macrozona!="RM") |> pull(macrozona))|>
  mutate(Dimension = factor(Dimension, levels = c("Macrozona Norte", "Macrozona Centro", "Macrozona Centro Sur", "Macrozona Sur", "Macrozona Austral")))|> 
  arrange(Dimension)|>
  knitr::kable("markdown", caption = "Comparación post-hoc, conglomerado-Macrozona")

Comparación post-hoc, conglomerado-Macrozona

Dimension	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
Macrozona Norte	4.79 (p<0.001)	-4.79 (p<0.001)
Macrozona Centro	-2.23 (p=0.259)	2.23 (p=0.259)
Macrozona Centro Sur	5.13 (p<0.001)	-5.13 (p<0.001)
Macrozona Sur	-7.19 (p<0.001)	7.19 (p<0.001)
Macrozona Austral	-0.19 (p=1.000)	0.19 (p=1.000)

pairwise_chisq_gof_test(dplyr::filter(tab_clus_macrozona_pam_om2_m,macrozona!="RM")[-1], p.adjust.method="holm")|>
  knitr::kable("markdown", caption="Dependencia categórica sol. 2 conglomerados (mensual), por pares de categorías en Macrozona (corrección Holm-Bonferroni)")
#0  ****
#Groups sharing a letter are not significantlt different (alpha = 0.05)

Dependencia categórica sol. 2 conglomerados (mensual), por pares de categorías en Macrozona (corrección Holm-Bonferroni)

n	group1	group2	statistic	p	df	p.adj	p.adj.signif
3682	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)	78.49767	0	4	0	****

2.1.3.f. Comparación covariables- Sexo

tab_clus_sexo_pam_om2_m<-  
  ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
    janitor::tabyl(glosa_sexo, clus_pam_om2 )

chisq_cramerv(tab_clus_sexo_pam_om2_m[,-1])
# $chisq_statistic
# [1] "205.87"
# 
# $chisq_df
# df 
#  1 
# 
# $chisq_p_value
# [1] "<0.001"
# 
# $cramers_v
# [1] "0.18"
#  
cat("Descartando valores negativos en sil width")
chisq_cramerv(with(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m), table(glosa_sexo, clus_pam_om2))
)#$chisq_p_value
#[1] "<0.001"

tab_clus_sexo_pam_om2_m|>
  janitor::adorn_percentages("col")|>
  janitor::adorn_rounding(digits = 3)|>
  dplyr::mutate(across(-glosa_sexo, ~sprintf("%d (%.1f%%)", tab_clus_sexo_pam_om2_m[[cur_column()]],.*100)))|> 
  (\(df) {
        if (interactive()) {mutate(df, across(-1, ~ gsub("%", "", gsub("\\.", ",", .))))|>  slice(2) |> rio::export("clipboard")}
        knitr::kable(df, "markdown", caption="Porcentajes por columna, conglomerado vs. sexo")
      })()

$chisq_statistic
[1] "205.87"

$chisq_df
df 
 1 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.18"

Descartando valores negativos en sil width$chisq_statistic
[1] "217.43"

$chisq_df
df 
 1 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.18"

Porcentajes por columna, conglomerado vs. sexo

glosa_sexo	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
HOMBRE	2466 (42.1%)	535 (69.6%)
MUJER	3391 (57.9%)	234 (30.4%)

2.1.3.g. Comparación covariables- Edad

dt_ing_calendar_month_t_desde_primera_adm_dedup|>
    dplyr::filter(month == 0)|>
    dplyr::inner_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run","clus_pam_om2")], by="run")|>
    dplyr::group_by(clus_pam_om2)|>
    dplyr::summarise(mean_edad = mean(min_edad_anos),
                     sd= sd(min_edad_anos),
                     p50= quantile(min_edad_anos,.5),
                     ci_lower = quantile(min_edad_anos, 0.25),
                     ci_upper = quantile(min_edad_anos, 0.75)) |> 
    dplyr::mutate(
    edad_ingreso = sprintf(
      "%.2f (%.2f ± %.2f); %.0f [%.0f, %.0f]",
      mean_edad, mean_edad, sd, p50, ci_lower, ci_upper
    )
  )|>
  dplyr::select(clus_pam_om2, edad_ingreso) |> 
    #dplyr::mutate_at(2:ncol(.),~scales::percent(.))|> 
    knitr::kable(caption="Descriptivos, edad minima de ingreso por conglomerado", digits=2, col.names=c("Conglomerado", "Promedio Desv. Estándar Mediana p025 p975"))

Descriptivos, edad minima de ingreso por conglomerado

Conglomerado	Promedio Desv. Estándar Mediana p025 p975
6623, Primer mes, TSM y Comorbilidad(2)	20.58 (20.58 ± 4.33); 20 [17, 24]
6612, Primer mes TUS(1)	22.38 (22.38 ± 4.32); 23 [18, 26]

dt_ing_calendar_month_t_desde_primera_adm_dedup|>
  dplyr::filter(month == 0)|>
  dplyr::inner_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run","clus_pam_om2")], by="run")|>
  dplyr::group_by(clus_pam_om2)|>
  dplyr::summarise(mean_edad = mean(min_edad_anos),
                   sd= sd(min_edad_anos),
                   ci_lower = mean(min_edad_anos) - qt(0.975, n()-1) * sd(min_edad_anos)/sqrt(n()),
                   ci_upper = mean(min_edad_anos) + qt(0.975, n()-1) * sd(min_edad_anos)/sqrt(n()))|> # Plot con ggplot2
ggplot(aes(x = clus_pam_om2, y = mean_edad)) +
  geom_point() +
  geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
  labs(title = NULL,
       x = "Conglomerado",
       y = "Edad promedio") +
  theme_minimal()+
  coord_flip()+
  theme(#axis.text.x = element_text(angle = 45, hjust = 1),
        panel.grid = element_blank())+
  theme(
    axis.text.y = element_text(size = 17, face = "bold"),#,margin = margin(l = 7)),           # Tamaño de las etiquetas de los grupos étnicos
    axis.text.x = element_text(size = 17, face = "bold"),           # Tamaño de las etiquetas del eje X
    axis.title.x = element_text(size = 16, face = "bold"),#,margin = margin(t = -15)),          # Tamaño del título del eje X
    axis.title.y = element_text(size = 16, face = "bold"),          # Tamaño del título del eje Y
    plot.title = NULL,  # Tamaño y estilo del título del gráfico
    legend.title = element_text(size = 17, face = "bold"),          # Tamaño del título de la leyenda
    legend.spacing.y = unit(1.5, "lines"),
    legend.box.spacing = unit(0.5, "lines"),      # Controla el espacio entre la leyenda y el gráfico
    legend.margin = margin(5, 5, 5, 5),  
    legend.key.height = unit(1, "cm"),  
    legend.text = element_text(size = 15, face = "bold")            # Tamaño del texto de la leyenda
  ) 
ggsave("_figs/edad_minima_por_cluster_pam_om2_m_25.png", dpi=600)

Edad promedio primer ingreso con intervalo de confianza por conglomerado

t.test(min_edad_anos ~ clus_pam_om2, 
             data = dt_ing_calendar_month_t_desde_primera_adm_dedup|>
               dplyr::filter(month == 0)|>
               dplyr::inner_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[,c("run","clus_pam_om2")], by="run"),var.equal = T)


dt_ing_calendar_month_t_desde_primera_adm_dedup|>
    dplyr::filter(month == 0)|>
    dplyr::inner_join(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens[, c("run", "clus_pam_om2")], by = "run")|>
   (\(df) {
        effsize::cohen.d(
     df$min_edad_anos,
     df$clus_pam_om2)
      })()
  
# d estimate: -0.4164028 (small)
# 95 percent confidence interval:
#      lower      upper 
# -0.4919251 -0.3408805 


cat("Descartando valores negativos en sil width")
dt_ing_calendar_month_t_desde_primera_adm_dedup|>
    dplyr::filter(month == 0)|>
    dplyr::inner_join(slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m)[,c("run","clus_pam_om2")], by="run")|>
   (\(df) {
    effsize::cohen.d(
     df$min_edad_anos,
     df$clus_pam_om2)
  })()

    Two Sample t-test

data:  min_edad_anos by clus_pam_om2
t = -10.856, df = 6624, p-value < 2.2e-16
alternative hypothesis: true difference in means between group 6623, Primer mes, TSM y Comorbilidad(2) and group 6612, Primer mes TUS(1) is not equal to 0
95 percent confidence interval:
 -2.125952 -1.475625
sample estimates:
mean in group 6623, Primer mes, TSM y Comorbilidad(2) 
                                             20.58153 
                mean in group 6612, Primer mes TUS(1) 
                                             22.38231 


Cohen's d

d estimate: -0.4164028 (small)
95 percent confidence interval:
     lower      upper 
-0.4919251 -0.3408805 

Descartando valores negativos en sil width
Cohen's d

d estimate: -0.4270268 (small)
95 percent confidence interval:
     lower      upper 
-0.5027064 -0.3513472 

2.1.3.h. Comparación covariables- Previsión

tab_clus_prev_pam_om2_m<-
ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|>
    janitor::tabyl(prev_benef_rec_post, clus_pam_om2)

tab_clus_prev_pam_om2_m|> 
     (\(df) {
     print(janitor::chisq.test(df))
     print(janitor::fisher.test(df, simulate.p.value = T, B = 1e5))
  })()

#X-squared = 15.684, df = 4, p-value = 0.003474
#p-value = 0.00189

chisq_cramerv(tab_clus_prev_pam_om2_m[-1])
# $chisq_statistic
# [1] "15.68"
# 
# $chisq_df
# df 
#  4 
# 
# $chisq_p_value
# [1] "0.0035"
# 
# $cramers_v
# [1] "0.05"
slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m)|>
    janitor::tabyl(prev_benef_rec_post, clus_pam_om2)|> 
     (\(df) {
     print(chisq_cramerv(df[-1]))
      print(janitor::chisq.test(df))
     print(janitor::fisher.test(df, simulate.p.value = T, B = 1e5))
  })()
# $chisq_p_value
# [1] "0.0024"
#p-value = 0.00147

tab_clus_prev_pam_om2_m|> 
  janitor::adorn_percentages("col")|> 
  janitor::adorn_rounding(digits = 3)|> 
  dplyr::mutate(across(-prev_benef_rec_post, ~sprintf("%d (%.1f%%)", tab_clus_prev_pam_om2_m[[cur_column()]],.*100)))|>
  #dplyr::mutate_at(2:ncol(.),~scales::percent(.))|> 
  knitr::kable(caption="Porcentajes por columna, conglomerado vs. Beneficios")

    Pearson's Chi-squared test

data:  df
X-squared = 15.684, df = 4, p-value = 0.003474


    Fisher's Exact Test for Count Data with simulated p-value (based on
    1e+05 replicates)

data:  df
p-value = 0.00225
alternative hypothesis: two.sided

$chisq_statistic
[1] "15.68"

$chisq_df
df 
 4 

$chisq_p_value
[1] "0.0035"

$cramers_v
[1] "0.05"

$chisq_statistic
[1] "16.54"

$chisq_df
df 
 4 

$chisq_p_value
[1] "0.0024"

$cramers_v
[1] "0.05"


    Pearson's Chi-squared test

data:  df
X-squared = 16.543, df = 4, p-value = 0.002371


    Fisher's Exact Test for Count Data with simulated p-value (based on
    1e+05 replicates)

data:  df
p-value = 0.00129
alternative hypothesis: two.sided

Porcentajes por columna, conglomerado vs. Beneficios

prev_benef_rec_post	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
FFAA	235 (4.0%)	14 (1.8%)
FONASA A	1554 (26.5%)	227 (29.5%)
FONASA BC	1902 (32.5%)	248 (32.2%)
FONASA D	673 (11.5%)	105 (13.7%)
ISAPRE	1493 (25.5%)	175 (22.8%)

2.1.3.i. Comparación covariables- Niv. Complejidad

tab_clus_compl_pam_om2_m<-
ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens|> 
  dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
  dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], 
                    by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first")|> 
  janitor::tabyl(nivel_de_complejidad, clus_pam_om2) 

chisq_cramerv(tab_clus_compl_pam_om2_m[-1])

janitor::chisq.test(tab_clus_compl_pam_om2_m)

cat("Sin siluetas negativas")
slice(ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens, -sil_neg_pam_om_clus2_m)|> 
  dplyr::mutate(estab_homo_base= as.numeric(estab_homo_base))|>
  dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], 
                    by = c("estab_homo_base" = "ESTAB_HOMO"), multiple = "first")|> 
  janitor::tabyl(nivel_de_complejidad, clus_pam_om2)|> 
  (\(df) {
          print(chisq_cramerv(df[-1]))
         print(janitor::chisq.test(df))
         print(janitor::fisher.test(df, simulate.p.value = T, B = 1e5))   
  })()
# $chisq_p_value
# [1] "<0.001"
# X-squared = 46.192, df = 4, p-value = 2.247e-09
# p-value = 1e-05

tab_clus_compl_pam_om2_m|> 
  janitor::adorn_percentages("col")|> 
  janitor::adorn_rounding(digits = 3)|> 
  dplyr::mutate(across(-nivel_de_complejidad, ~sprintf("%d (%.1f%%)", tab_clus_compl_pam_om2_m[[cur_column()]],.*100)))|>
  #dplyr::mutate_at(2:ncol(.), ~scales::percent(as.numeric(.), accuracy=.1))|> 
  knitr::kable(caption="Tabla de contingencia, Niv. de complejidad (proporción por columna)")

$chisq_statistic
[1] "45.00"

$chisq_df
df 
 4 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.08"


    Pearson's Chi-squared test

data:  tab_clus_compl_pam_om2_m
X-squared = 45.003, df = 4, p-value = 3.97e-09

Sin siluetas negativas$chisq_statistic
[1] "46.19"

$chisq_df
df 
 4 

$chisq_p_value
[1] "<0.001"

$cramers_v
[1] "0.08"


    Pearson's Chi-squared test

data:  df
X-squared = 46.192, df = 4, p-value = 2.247e-09


    Fisher's Exact Test for Count Data with simulated p-value (based on
    1e+05 replicates)

data:  df
p-value = 1e-05
alternative hypothesis: two.sided

Tabla de contingencia, Niv. de complejidad (proporción por columna)

nivel_de_complejidad	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
Alta Complejidad	3547 (60.6%)	520 (67.6%)
Baja Complejidad	940 (16.0%)	147 (19.1%)
Mediana Complejidad	1160 (19.8%)	91 (11.8%)
Pendiente	167 (2.9%)	5 (0.7%)
Sin dato	43 (0.7%)	6 (0.8%)

chisq.posthoc.test(tab_clus_compl_pam_om2_m[-1])|> 
 (\(df) dplyr::mutate(df, across(3:ncol(df), ~round(as.numeric(gsub("\\*","",.x)),3))))()|>
  # knitr::kable("html", caption="Comparación post-hoc, conglomerado-región") |>
  (\(df) dplyr::group_by(df, Dimension))()|>
  (\(df) dplyr::summarise(df, across(2:(ncol(df) - 1), ~paste0(first(sprintf("%1.2f", .x)), " (p=", last(sprintf("%1.3f", .x)), ")"))))()|>
  (\(df) dplyr::mutate(df, across(2:length(names(df)), ~gsub("p=0.000)", "p<0.001)",.x))))()|>
  ungroup()|> 
  dplyr::mutate(Dimension = dplyr::pull(tab_clus_compl_pam_om2_m[1]))|>
  mutate(Dimension = factor(Dimension, levels = c("Alta Complejidad", "Mediana Complejidad", "Baja Complejidad", "Pendiente", "Sin dato")))|> 
  arrange(Dimension)|>
  (\(df) {
    if (interactive()) {mutate(df, across(-Dimension, ~ gsub("\\(|\\)", "", gsub("\\.", ",", .))))|> rio::export("clipboard")}
    knitr::kable(df, "markdown", caption="Comparación post-hoc, conglomerado-Niv. complejidad")
  })() 

Comparación post-hoc, conglomerado-Niv. complejidad

Dimension	6623, Primer mes, TSM y Comorbilidad(2)	6612, Primer mes TUS(1)
Alta Complejidad	-3.78 (p=0.002)	3.78 (p=0.002)
Mediana Complejidad	5.31 (p<0.001)	-5.31 (p<0.001)
Baja Complejidad	-2.16 (p=0.308)	2.16 (p=0.308)
Pendiente	3.61 (p=0.003)	-3.61 (p=0.003)
Sin dato	-0.14 (p=1.000)	0.14 (p=1.000)

2.1.4. Compilación comparación covariables

#dput(attr(t(tab_clus_compl_pam_om7_q),"dimnames")[[1]])

# Definir los datos correctamente
data_pam_om2_m <-cbind.data.frame(
  Grupo= c("6035, Un trimestre, TSM(2)", "6025, Un trimestre, TUS(1)"), 
      PPOO_bin            = c(NA, NA), 
      PPOO_sinautoid          = c(NA, NA), 
      PPOO_conautoid          = c(NA, NA), 
      Mortalidad              = c(NA, "+"), 
      RM                      = c(NA, "-"), 
      `Macrozona-Austral`       = c(NA, NA), 
      `Macrozona-Centro`        = c(NA, "-"), 
      `Macrozona-Centro Sur`    = c(NA, NA), 
      `Macrozona-Norte`        = c(NA, "-"), 
      `Macrozona-Sur`         = c(NA, "+"), 
      Sexo_mujeres             = c(NA, "-"), 
      `Edad ingreso`         = c(NA, "+"), 
      `Previsión-FFAA`         = c(NA, NA), 
      `Previsión-FONASA A`     = c(NA, NA), 
      `Previsión-FONASA BC`     = c(NA, NA), 
      `Previsión-FONASA D`     = c(NA, NA), 
      `Previsión-ISAPRE`     = c(NA, NA), 
      `NivComp-Baja`     = c(NA, NA), 
      `NivComp-Media`     = c(NA, "+"), 
      `NivComp-Alta`     = c(NA, "-"))
#
# Asegurar que los nombres de las columnas sean válidos y no haya espacios en blanco
# Derretir el dataframe para que sea adecuado para ggplot2
data_melt_pam_om2_m <- reshape2::melt(data_pam_om2_m, id.vars = 'Grupo', variable.name = 'Variable', value.name = 'Asociación')

# Reemplazar los NA por un valor vacío
data_melt_pam_om2_m$Asociación[is.na(data_melt_pam_om2_m$Asociación)] <- "\n"

# Crear el gráfico con ggplot
data_melt_pam_om2_m|> 
  dplyr::filter(Grupo=="6025, Un trimestre, TUS(1)")|>
  dplyr::mutate(Variable = gsub("_", " ", Variable))|> 
ggplot(aes(x = Variable, y = Grupo, fill = Asociación)) +
  geom_tile(color = "white", size = 0.8) +
  scale_fill_manual(values = c("+" = "#556B2F", "-" = "#E2725B", "\n" = "white")) +
  labs(title =NULL, x = "Variables", y = "Conglomerado") +
  theme_minimal() +
  theme(#axis.text.x = element_text(angle = 45, hjust = 1),
        panel.grid = element_blank())+
  theme(
    axis.text.y = element_text(size = 17, face = "bold"),#,margin = margin(l = 7)),           # Tamaño de las etiquetas de los grupos étnicos
    axis.text.x = element_text(size = 17, face = "bold"),           # Tamaño de las etiquetas del eje X
    axis.title.x = element_text(size = 16, face = "bold"),#,margin = margin(t = -15)),          # Tamaño del título del eje X
    axis.title.y = element_text(size = 16, face = "bold"),          # Tamaño del título del eje Y
    plot.title = NULL,  # Tamaño y estilo del título del gráfico
    legend.title = element_text(size = 17, face = "bold"),          # Tamaño del título de la leyenda
    legend.spacing.y = unit(1.5, "lines"),
    legend.box.spacing = unit(0.5, "lines"),      # Controla el espacio entre la leyenda y el gráfico
    legend.margin = margin(5, 5, 5, 5),  
    legend.key.height = unit(1, "cm"),  
    legend.text = element_text(size = 15, face = "bold")            # Tamaño del texto de la leyenda
  ) +
  coord_flip()

ggsave("_figs/asociaciones_pam_om2_m_25.png", width=8.8*.8, height=5*.8, dpi=1000)

Información de la sesión

cat(paste0("R library: ", Sys.getenv("R_LIBS_USER")))

R library: C:\Users\andre\AppData\Local/R/win-library/4.4

cat(paste0("Date: ",withr::with_locale(new = c('LC_TIME' = 'C'), code =Sys.time())))

Date: 2025-05-13 22:16:32.702111

cat(paste0("Editor context: ", getwd()))

Editor context: H:/My Drive/PERSONAL ANDRES/UCH_salud_publica/asignaturas/un_inv_II

cat("quarto version: "); system("quarto --version") 

quarto version: 

[1] 0

save.image("avance250117_3_25.RData")

sesion_info <- devtools::session_info()

Warning in system2("quarto", "-V", stdout = TRUE, env = paste0("TMPDIR=", : el
comando ejecutado '"quarto"
TMPDIR=C:/Users/andre/AppData/Local/Temp/RtmpwVf96u/file148c490be4a -V' tiene
el estatus 1

dplyr::select(
  tibble::as_tibble(sesion_info$packages),
  c(package, loadedversion, source)
)|> 
 knitr::kable(caption = "R packages", format = "html",
      col.names = c("Row number", "Package", "Version"),
    row.names = FALSE,
      align = c("c", "l", "r"))|> 
  kableExtra::kable_styling(bootstrap_options = c("striped", "hover"),font_size = 12)|> 
  kableExtra::scroll_box(width = "100%", height = "375px")  

R packages

Row number	Package	Version
abind	1.4-8	CRAN (R 4.4.1)
backports	1.5.0	CRAN (R 4.4.1)
biostat3	0.2.2	CRAN (R 4.4.3)
boot	1.3-31	CRAN (R 4.4.1)
broom	1.0.7	CRAN (R 4.4.1)
cachem	1.1.0	CRAN (R 4.4.1)
car	3.1-3	CRAN (R 4.4.1)
carData	3.0-5	CRAN (R 4.4.1)
chisq.posthoc.test	0.1.3	Github (ebbertd/chisq.posthoc.test@186d2ca6bbdba9fc19601aff4696ae1b85e7e0b0)
cli	3.6.4	CRAN (R 4.4.1)
cluster	2.1.8.1	CRAN (R 4.4.1)
coda	0.19-4.1	CRAN (R 4.4.1)
codetools	0.2-20	CRAN (R 4.4.1)
colorspace	2.1-1	CRAN (R 4.4.1)
cowplot	1.1.3	CRAN (R 4.4.1)
crayon	1.5.3	CRAN (R 4.4.1)
curl	6.2.1	CRAN (R 4.4.1)
data.table	1.17.0	CRAN (R 4.4.1)
devtools	2.4.5	CRAN (R 4.4.1)
DiagrammeR	1.0.11	CRAN (R 4.4.1)
DiagrammeRsvg	0.1	CRAN (R 4.4.1)
digest	0.6.37	CRAN (R 4.4.1)
doFuture	1.0.1	CRAN (R 4.4.1)
doParallel	1.0.17	CRAN (R 4.4.1)
dplyr	1.1.4	CRAN (R 4.4.3)
effsize	0.8.1	CRAN (R 4.4.1)
ellipsis	0.3.2	CRAN (R 4.4.1)
emmeans	1.11.0	CRAN (R 4.4.3)
epitools	0.5-10.1	CRAN (R 4.4.0)
estimability	1.5.1	CRAN (R 4.4.1)
evaluate	1.0.3	CRAN (R 4.4.1)
expsmooth	2.3	CRAN (R 4.4.3)
factoextra	1.0.7	CRAN (R 4.4.1)
farver	2.1.2	CRAN (R 4.4.1)
fastcluster	1.2.6	CRAN (R 4.4.1)
fastmap	1.2.0	CRAN (R 4.4.1)
fma	2.5	CRAN (R 4.4.3)
forcats	1.0.0	CRAN (R 4.4.1)
foreach	1.5.2	CRAN (R 4.4.1)
forecast	8.24.0	CRAN (R 4.4.3)
Formula	1.2-5	CRAN (R 4.4.1)
fpp2	2.5	CRAN (R 4.4.3)
fracdiff	1.5-3	CRAN (R 4.4.1)
fs	1.6.5	CRAN (R 4.4.1)
future	1.34.0	CRAN (R 4.4.1)
future.apply	1.11.3	CRAN (R 4.4.1)
generics	0.1.3	CRAN (R 4.4.1)
ggh4x	0.3.0	CRAN (R 4.4.1)
ggplot2	3.5.2	CRAN (R 4.4.3)
ggpubr	0.6.0	CRAN (R 4.4.1)
ggrepel	0.9.6	CRAN (R 4.4.1)
ggseqplot	0.8.5	CRAN (R 4.4.1)
ggsignif	0.6.4	CRAN (R 4.4.1)
globals	0.16.3	CRAN (R 4.4.1)
glue	1.8.0	CRAN (R 4.4.1)
gridExtra	2.3	CRAN (R 4.4.1)
gtable	0.3.6	CRAN (R 4.4.1)
gtsummary	2.2.0	CRAN (R 4.4.3)
haven	2.5.4	CRAN (R 4.4.1)
hms	1.1.3	CRAN (R 4.4.1)
htmltools	0.5.8.1	CRAN (R 4.4.1)
htmlwidgets	1.6.4	CRAN (R 4.4.1)
httpuv	1.6.15	CRAN (R 4.4.1)
iterators	1.0.14	CRAN (R 4.4.1)
janitor	2.2.1	CRAN (R 4.4.1)
job	0.3.1	CRAN (R 4.4.1)
jsonlite	1.9.1	CRAN (R 4.4.1)
kableExtra	1.4.0	CRAN (R 4.4.1)
km.ci	0.5-6	CRAN (R 4.4.3)
KMsurv	0.1-5	CRAN (R 4.4.0)
knitr	1.49	CRAN (R 4.4.1)
labeling	0.4.3	CRAN (R 4.4.1)
later	1.4.1	CRAN (R 4.4.1)
lattice	0.22-6	CRAN (R 4.4.1)
lifecycle	1.0.4	CRAN (R 4.4.1)
listenv	0.9.1	CRAN (R 4.4.1)
lmtest	0.9-40	CRAN (R 4.4.1)
lubridate	1.9.4	CRAN (R 4.4.1)
magrittr	2.0.3	CRAN (R 4.4.1)
MASS	7.3-60.2	CRAN (R 4.4.1)
Matrix	1.7-0	CRAN (R 4.4.1)
memoise	2.0.1	CRAN (R 4.4.1)
mgcv	1.9-1	CRAN (R 4.4.1)
mime	0.12	CRAN (R 4.4.1)
miniUI	0.1.1.1	CRAN (R 4.4.1)
mnormt	2.1.1	CRAN (R 4.4.1)
multcomp	1.4-28	CRAN (R 4.4.1)
munsell	0.5.1	CRAN (R 4.4.1)
mvtnorm	1.3-3	CRAN (R 4.4.1)
NbClust	3.0.1	CRAN (R 4.4.1)
nlme	3.1-164	CRAN (R 4.4.1)
nnet	7.3-19	CRAN (R 4.4.1)
pacman	0.5.1	CRAN (R 4.4.3)
parallelly	1.42.0	CRAN (R 4.4.1)
patchwork	1.3.0	CRAN (R 4.4.1)
permute	0.9-7	CRAN (R 4.4.1)
pillar	1.10.1	CRAN (R 4.4.1)
pkgbuild	1.4.6	CRAN (R 4.4.1)
pkgconfig	2.0.3	CRAN (R 4.4.1)
pkgload	1.4.0	CRAN (R 4.4.1)
plyr	1.8.9	CRAN (R 4.4.1)
profvis	0.4.0	CRAN (R 4.4.1)
progressr	0.15.1	CRAN (R 4.4.1)
promises	1.3.2	CRAN (R 4.4.1)
psych	2.4.12	CRAN (R 4.4.1)
purrr	1.0.4	CRAN (R 4.4.1)
quadprog	1.5-8	CRAN (R 4.4.1)
quantmod	0.4.26	CRAN (R 4.4.1)
R.methodsS3	1.8.2	CRAN (R 4.4.1)
R.oo	1.27.0	CRAN (R 4.4.1)
R.utils	2.13.0	CRAN (R 4.4.1)
R6	2.6.1	CRAN (R 4.4.1)
ragg	1.3.3	CRAN (R 4.4.1)
rbibutils	2.3	CRAN (R 4.4.1)
RColorBrewer	1.1-3	CRAN (R 4.4.0)
Rcpp	1.0.14	CRAN (R 4.4.1)
Rdpack	2.6.2	CRAN (R 4.4.1)
readr	2.1.5	CRAN (R 4.4.3)
remotes	2.5.0	CRAN (R 4.4.1)
reshape2	1.4.4	CRAN (R 4.4.1)
rio	1.2.3	CRAN (R 4.4.1)
rlang	1.1.5	CRAN (R 4.4.1)
rmarkdown	2.29	CRAN (R 4.4.1)
rstatix	0.7.2	CRAN (R 4.4.1)
rstudioapi	0.17.1	CRAN (R 4.4.1)
rsvg	2.6.2	CRAN (R 4.4.3)
sandwich	3.1-1	CRAN (R 4.4.1)
scales	1.3.0	CRAN (R 4.4.1)
sessioninfo	1.2.3	CRAN (R 4.4.1)
shiny	1.10.0	CRAN (R 4.4.1)
snakecase	0.11.1	CRAN (R 4.4.1)
stargazer	5.2.3	CRAN (R 4.4.1)
stringi	1.8.4	CRAN (R 4.4.1)
stringr	1.5.1	CRAN (R 4.4.3)
survival	3.6-4	CRAN (R 4.4.1)
survminer	0.5.0	CRAN (R 4.4.3)
survMisc	0.5.6	CRAN (R 4.4.3)
svglite	2.1.3	CRAN (R 4.4.1)
systemfonts	1.2.1	CRAN (R 4.4.1)
textshaping	1.0.0	CRAN (R 4.4.1)
TH.data	1.1-3	CRAN (R 4.4.1)
tibble	3.2.1	CRAN (R 4.4.1)
tidyr	1.3.1	CRAN (R 4.4.1)
tidyselect	1.2.1	CRAN (R 4.4.1)
tidyverse	2.0.0	CRAN (R 4.4.3)
timechange	0.3.0	CRAN (R 4.4.1)
timeDate	4041.110	CRAN (R 4.4.1)
Tmisc	1.0.1	CRAN (R 4.4.1)
TraMineR	2.2-11	CRAN (R 4.4.1)
TraMineRextras	0.6.8	CRAN (R 4.4.1)
tseries	0.10-58	CRAN (R 4.4.1)
TTR	0.24.4	CRAN (R 4.4.1)
tzdb	0.4.0	CRAN (R 4.4.1)
urca	1.3-4	CRAN (R 4.4.1)
urlchecker	1.0.1	CRAN (R 4.4.1)
usethis	3.1.0	CRAN (R 4.4.1)
utf8	1.2.4	CRAN (R 4.4.1)
V8	6.0.1	CRAN (R 4.4.1)
vctrs	0.6.5	CRAN (R 4.4.1)
vegan	2.6-10	CRAN (R 4.4.1)
vegclust	2.0.2	CRAN (R 4.4.1)
viridisLite	0.4.2	CRAN (R 4.4.1)
visNetwork	2.1.2	CRAN (R 4.4.1)
WeightedCluster	1.8-1	CRAN (R 4.4.1)
withr	3.0.2	CRAN (R 4.4.1)
writexl	1.5.1	CRAN (R 4.4.1)
xfun	0.51	CRAN (R 4.4.1)
xml2	1.3.7	CRAN (R 4.4.1)
xtable	1.8-4	CRAN (R 4.4.1)
xts	0.14.1	CRAN (R 4.4.1)
yaml	2.3.10	CRAN (R 4.4.1)
zoo	1.8-13	CRAN (R 4.4.1)

reticulate::py_list_packages()|> 
 knitr::kable(caption = "Python packages", format = "html",
      col.names = c("Package", "Version", "Requirement"),
    row.names = FALSE,
      align = c("c", "l", "r", "r"))|> 
  kableExtra::kable_styling(bootstrap_options = c("striped", "hover"),font_size = 12)|> 
  kableExtra::scroll_box(width = "100%", height = "375px")  

Error in path.expand(path): argumento 'path' inválido